void plots( char *store)
{
 gROOT->Reset();
 gROOT->SetStyle("Plain");
 gStyle->SetOptStat(1);
 gStyle->SetTitleFont(72);
 TLatex *l = new TLatex();
 l->SetTextFont(72);
 TLine *line = new TLine();

 char  title[500];
 char  txtdir[512];
 
 // char *store = "files.list";

 TString outdir( gSystem->Getenv( "TOFRESULTDIR" ) );
 char* coutdir=char( outdir );
 sprintf(txtdir,"%s/files",coutdir);

 cout << "Input filename is " << store << endl;
 const int NumBars = 216;

 double rese[NumBars], reseErr[NumBars], resw[NumBars], reswErr[NumBars], 
        slope_rese[NumBars], slope_reseErr[NumBars], slope_resw[NumBars], slope_reswErr[NumBars];
 double speed[NumBars], speedErr[NumBars], speed_off[NumBars], speed_offErr[NumBars];
 double spdtw[NumBars], spdtwErr[NumBars], spdtw_off[NumBars], spdtw_offErr[NumBars];
 double twalke[NumBars], twalkeErr[NumBars], twalke_off[NumBars], twalke_offErr[NumBars];
 double twalkw[NumBars], twalkwErr[NumBars], twalkw_off[NumBars], twalkw_offErr[NumBars];
 double offsete[NumBars], offseteErr[NumBars], offsete_off[NumBars], offsete_offErr[NumBars];
 double offsetw[NumBars], offsetwErr[NumBars], offsetw_off[NumBars], offsetw_offErr[NumBars];
 double zrese[6][NumBars], zresw[6][NumBars];
 double q_mpv_e[NumBars], lambda_e[NumBars], a_e[NumBars], sigma_e[NumBars];  
 double q_mpv_w[NumBars], lambda_w[NumBars], a_w[NumBars], sigma_w[NumBars];  
 double dummy[27][NumBars];  
 
 TH1F *h1f_res         = new TH1F("h1f_res","PMT Resolutions", 100, 0, 300);
 TH1F *h1f_res_err     = new TH1F("h1f_res_err","PMT Resolutions errors", 100, 0, 60);
 TH1F *h1f_slope_res   = new TH1F("h1f_slope_res","PMT Slope Resolutions", 100, 0, 4);
 TH1F *h1f_slope_res_err = new TH1F("h1f_slope_res_err","PMT Slope Resolutions errors", 100, 0, 0.4);
 TH1F *h1f_speed       = new TH1F("h1f_speed","Speed of light before twalk", 100, 12, 16);
 TH1F *h1f_speed_err   = new TH1F("h1f_speed_err","Speed of light before twalk error", 100, 0, 0.05);
 TH1F *h1f_spdtw       = new TH1F("h1f_spdtw","Speed of light after twalk", 100, 13, 17);
 TH1F *h1f_spdtw_err   = new TH1F("h1f_spdtw_err","Speed of light after twalk error", 100, 0, 0.05);
 TH1F *h1f_twalke      = new TH1F("h1f_twalke","Time-walk (East side)", 100, 0, 100);
 TH1F *h1f_twalkw      = new TH1F("h1f_twalkw","Time-walk (West side)", 100, 0, 100);
 TH1F *h1f_twalke_err  = new TH1F("h1f_twalke_err","Time-walk Error (East side)", 100, 0, 6.2);
 TH1F *h1f_twalkw_err  = new TH1F("h1f_twalkw_err","Time-walk Error (West side)", 100, 0, 5.2);
 TH1F *h1f_offsete     = new TH1F("h1f_offsete","Offsets (East side)", 100, 52, 70);
 TH1F *h1f_offsetw     = new TH1F("h1f_offsetw","Offsets (West side)", 100, 52, 70);
 TH1F *h1f_offsete_err = new TH1F("h1f_offsete_err","Offset err (East side)", 100, 0, 3); 
 TH1F *h1f_offsetw_err = new TH1F("h1f_offsetw_err","Offset err (West side)", 100, 0, 3);
 TH1F *h1f_zres_0      = new TH1F("h1f_zres_0","Res ~ct", 100, -3, 3);
 TH1F *h1f_zres_1      = new TH1F("h1f_zres_1","Res ~z", 100, -0.008, 0.008);
 TH1F *h1f_zres_2      = new TH1F("h1f_zres_2","Res ~z*z", 100, -0.00002, 0.00002);
 TH1F *h1f_zres_3      = new TH1F("h1f_zres_3","Res ~1/sqrt(Q)", 100, -80, 80);
 TH1F *h1f_zres_4      = new TH1F("h1f_zres_4","Res ~z/sqrt(Q)", 100, -0.2, 0.2);
 TH1F *h1f_zres_5      = new TH1F("h1f_zres_5","Res ~1/Q", 100, -600, 600);
 TH1F *h1f_qmpv        = new TH1F("h1f_qmpv","Q MPV", 100, 0, 3000);
 TH1F *h1f_lambda      = new TH1F("h1f_lambda","Lambda", 100, 0, 400);
 TH1F *h1f_a           = new TH1F("h1f_a","a", 100, 0, 0.00003);
 TH1F *h1f_sigma       = new TH1F("h1f_sigma","sigma Q MPV", 100, 0, 400);
 TH2F *h2f_resol_qmpv  = new TH2F("h2f_resol_qmpv","Resol vs Q MPV", 100, 0, 3000, 100, 0, 300);


 // TPostScript myps(title,111);

//--------------------------------------------------------------------------
// Speed of light before time walk
//--------------------------------------------------------------------------
// myps.NewPage();

 ifstream speed_f;
 sprintf(title,"%s/%s_speed-calib.txt",txtdir,store);
 speed_f.open(title);

 for (int i=0; i<NumBars; i++) {
    if (!speed_f.good()) break;
 	speed_f >> dummy[0][i]
 		>> speed[i]
 		>> speedErr[i]
 		>> speed_off[i]
 		>> speed_offErr[i]
        	>> dummy[1][i]
        	>> dummy[2][i]
        	>> dummy[3][i]
        	>> dummy[4][i];
 }
 speed_f.close();      

 for (int i=0; i<NumBars; i++) {
    if (speed[i]<12.1 || speed[i]>15.5){ 
       cout << "Bar " << i << ": speed = " << speed[i] 
            << " +- " << speedErr[i] << endl;
    }else{
    h1f_speed->Fill(speed[i]); 
    h1f_speed_err->Fill(speedErr[i]);
    }
    if (speedErr[i]<0.01 || speedErr[i]>0.04) 
       cout << "Bar " << i << ": speed err = " << speedErr[i] << endl;
 }

//--------------------------------------------------------------------------
// Speed of light after time walk
//--------------------------------------------------------------------------
// myps.NewPage();

 ifstream spdtw_f;
 sprintf(title,"%s/%s_spdtw-calib.txt",txtdir,store);
 spdtw_f.open(title);

 for (int i=0; i<NumBars; i++) {
    if (!spdtw_f.good()) break;
 	spdtw_f >> dummy[0][i]
 		>> spdtw[i]
 		>> spdtwErr[i]
 		>> spdtw_off[i]
 		>> spdtw_offErr[i]
        	>> dummy[1][i]
        	>> dummy[2][i]
        	>> dummy[3][i]
        	>> dummy[4][i];
 }
 spdtw_f.close();      

 for (int i=0; i<NumBars; i++) {
    if (spdtw[i]<14 || spdtw[i]>16.5){ 
       cout << "Bar " << i << ": spdtw = " << spdtw[i] 
            << " +- " << spdtwErr[i] << endl;
    }else{
    h1f_spdtw->Fill(spdtw[i]); 
    h1f_spdtw_err->Fill(spdtwErr[i]);
    }
    if (spdtwErr[i]<=0 || spdtwErr[i]>0.2) 
       cout << "Bar " << i << ": spdtw err = " << spdtwErr[i] << endl;
 }

 TCanvas *c2 = new TCanvas("c2","c2",800,1000);
 c2->SetFillColor(10);
 c2->SetBottomMargin(0.2);
 c2->SetLeftMargin(0.2);
 c2->Divide(2,2);

 c2->cd(1);
 Draw(h1f_speed,4);

 c2->cd(2);
 Draw(h1f_speed_err,4);

 c2->cd(3);
 Draw(h1f_spdtw,4);

 c2->cd(4);
 Draw(h1f_spdtw_err,4);

 // c2->Update();
 
 char speedt[100];
 sprintf(speedt,"%s_speed.eps",store);
 c2->Print(speedt);
//--------------------------------------------------------------------------
// Time-walk east and west and errors
//--------------------------------------------------------------------------
// myps.NewPage();

 ifstream twalke_f;
 sprintf(title,"%s/%s_timew-calib_e.txt",txtdir,store);
 twalke_f.open(title);

 for (int i=0; i<NumBars; i++) {
    if (!twalke_f.good()) break;
    twalke_f >> dummy[0][i]
 		  >> twalke[i]
 		  >> twalkeErr[i];
 }

 twalke_f.close();        

 for (int i=0; i<NumBars; i++) {
      if(twalke[i]<40 || twalke[i] > 100){
        cout << "PMT " << i << " (east): Twalk = " << twalke[i] 
	<< " +- " << twalkeErr[i] << endl;
      }else{
      h1f_twalke->Fill(twalke[i]); 
      h1f_twalke_err->Fill(twalkeErr[i]);   
      } 
      if(twalkeErr[i]<0 || twalkeErr[i] > 4)
        cout << "PMT " << i << " (east): Twalk err = " << twalkeErr[i] << endl;
 }  

 ifstream twalkw_f;
 sprintf(title,"%s/%s_timew-calib_w.txt",txtdir,store);
 twalkw_f.open(title);

 for (int i=0; i<NumBars; i++) {
    if (!twalkw_f.good()) break;
    twalkw_f >> dummy[0][i]
 		  >> twalkw[i]
 		  >> twalkwErr[i];
 }

 twalkw_f.close();        

 for (int i=0; i<NumBars; i++) {
      if(twalkw[i]<40 || twalkw[i] > 100){
        cout << "PMT " << i << " (west): Twalk = " << twalkw[i] 
	     << " +- " << twalkwErr[i] << endl;
      }else{
      h1f_twalkw->Fill(twalkw[i]);
      h1f_twalkw_err->Fill(twalkwErr[i]);
      }  
      if(twalkwErr[i]<=0 || twalkwErr[i] > 4)
        cout << "PMT " << i << " (west): Twalk err = " << twalkwErr[i] << endl;
 }  

 TCanvas *c3 = new TCanvas("c3","c3",800,1000);
 c3->SetFillColor(10);
 c3->SetBottomMargin(0.2);
 c3->SetLeftMargin(0.2);
 c3->Divide(2,2);

 c3->cd(1);
 Draw(h1f_twalke,4);

 c3->cd(2);
 Draw(h1f_twalkw,4);

 c3->cd(3);
 Draw(h1f_twalke_err,4);

 c3->cd(4);
 Draw(h1f_twalkw_err,4);

 //c3->Update();
 char walkt[100];
 sprintf(walkt,"%s_walk.eps",store);
 c3->Print(walkt);

//--------------------------------------------------------------------------
// Offsets east and west and errors
//--------------------------------------------------------------------------
// myps.NewPage();

 ifstream offsete_f;
 sprintf(title,"%s/%s_T0-calib_e.txt",txtdir,store);
 offsete_f.open(title);

 for (int i=0; i<NumBars; i++) {
    if (!offsete_f.good()) break;
    offsete_f >> dummy[0][i]
        	  >> dummy[1][i]
        	  >> dummy[2][i]
 		  >> offsete[i]
                  >> dummy[3][i];
 }

 offsete_f.close();        

 for (int i=0; i<NumBars; i++) {
      if(offsete[i]<45 || offsete[i]>65){
        cout << "PMT " << i << " (east): Offset = " << offsete[i] << endl;
      }else{
      h1f_offsete->Fill(offsete[i]);  
      }
 }  

 ifstream offsetw_f;
 sprintf(title,"%s/%s_T0-calib_w.txt",txtdir,store);
 offsetw_f.open(title);

 for (int i=0; i<NumBars; i++) {
    if (!offsetw_f.good()) break;
    offsetw_f >> dummy[0][i]
 		  >> offsetw[i]
        	  >> dummy[1][i]
        	  >> dummy[2][i]
        	  >> dummy[3][i];
 }

 offsetw_f.close();        

 for (int i=0; i<NumBars; i++) {
      if(offsetw[i]<45 || offsetw[i]>65){
        cout << "PMT " << i << " (west): Offset = " << offsetw[i] << endl;
      }else{
        h1f_offsetw->Fill(offsetw[i]);  
      }
 }  

//--------------------------------------------------------------------------
// Resolution
//--------------------------------------------------------------------------
// myps.NewPage();

 ifstream rese_f;
 sprintf(title,"%s/%s_PmtRes-calib_e.txt",txtdir,store);
 rese_f.open(title);

 for (int i=0; i<NumBars; i++) {
    if (!rese_f.good()) break;
 	rese_f  >> dummy[0][i]
 		>> rese[i]
 		>> reseErr[i]
		>> slope_rese[i]
        	>> slope_reseErr[i];
 }
 rese_f.close();      

 ifstream resw_f;
 sprintf(title,"%s/%s_PmtRes-calib_w.txt",txtdir,store);
 resw_f.open(title);

 for (int i=0; i<NumBars; i++) {
    if (!resw_f.good()) break;
 	resw_f  >> dummy[0][i]
 		>> resw[i]
 		>> reswErr[i]
		>> slope_resw[i]
        	>> slope_reswErr[i];
 }
 resw_f.close();      

 for (int i=0; i<NumBars; i++) {
    if (rese[i]*1000<40 || rese[i]*1000>220){ 
       cout << "PMT " << i << " (east): res = " << rese[i]*1000 << endl;
    }
       h1f_res->Fill(rese[i]*1000); 
       h1f_res_err->Fill(reseErr[i]*1000);

    if (resw[i]*1000<40 || resw[i]*1000>220){ 
       cout << "PMT " << i << " (west): res = " << resw[i]*1000 << endl;
    }
       h1f_res->Fill(resw[i]*1000); 
       h1f_res_err->Fill(reswErr[i]*1000);


    if (reseErr[i]*1000<=0 || reseErr[i]*1000>20) 
       cout << "PMT " << i << " (east): res err = " << reseErr[i]*1000 << endl;
    if (reswErr[i]*1000<=0 || reswErr[i]*1000>20) 
       cout << "PMT " << i << " (west): res err = " << reswErr[i]*1000 << endl;

    if (slope_rese[i]*1000<=0 || slope_rese[i]*1000>2){ 
       cout << "PMT " << i << " (east): slope res = " << slope_rese[i]*1000 << endl;
    }
       h1f_slope_res->Fill(slope_rese[i]*1000); 
       h1f_slope_res_err->Fill(slope_reseErr[i]*1000);

    if (slope_resw[i]*1000<=0 || slope_resw[i]*1000>2){ 
       cout << "PMT " << i << " (west): slope res = " << slope_resw[i]*1000 << endl;
    }
       h1f_slope_res->Fill(slope_resw[i]*1000); 
       h1f_slope_res_err->Fill(slope_reswErr[i]*1000);
 

    if (slope_reseErr[i]*1000<=0 || slope_reseErr[i]*1000>0.2) 
       cout << "PMT " << i << " (east): slope res err = " << slope_reseErr[i]*1000 << endl;
    if (slope_reswErr[i]*1000<=0 || slope_reswErr[i]*1000>0.2) 
       cout << "PMT " << i << " (west): slope res err = " << slope_reswErr[i]*1000 << endl;
 }
/*
 for (int i=0; i<NumBars; i++) {
    if (rese[i]*1000<40 || rese[i]*1000>220){ 
       cout << "PMT " << i << " (east): res = " << rese[i]*1000 << endl;
    }else{
       h1f_res->Fill(rese[i]*1000); 
       h1f_res_err->Fill(reseErr[i]*1000);
    }
    if (resw[i]*1000<40 || resw[i]*1000>220){ 
       cout << "PMT " << i << " (west): res = " << resw[i]*1000 << endl;
    }else{
       h1f_res->Fill(resw[i]*1000); 
       h1f_res_err->Fill(reswErr[i]*1000);
    }
    if (reseErr[i]*1000<=0 || reseErr[i]*1000>20) 
       cout << "PMT " << i << " (east): res err = " << reseErr[i]*1000 << endl;
    if (reswErr[i]*1000<=0 || reswErr[i]*1000>20) 
       cout << "PMT " << i << " (west): res err = " << reswErr[i]*1000 << endl;
 }
*/
 TCanvas *c0 = new TCanvas("c0","c0",800,1000);
 c0->SetFillColor(10);
 c0->SetBottomMargin(0.2);
 c0->SetLeftMargin(0.2);
 c0->Divide(2,3);
 
 c0->cd(1);
 Draw(h1f_offsete,4);

 c0->cd(2);
 Draw(h1f_offsetw,4);

 c0->cd(3);
 Draw(h1f_res,4);
 
 c0->cd(4);
 Draw(h1f_res_err,4);

 c0->cd(5);
 Draw(h1f_slope_res,4);
 
 c0->cd(6);
 Draw(h1f_slope_res_err,4);
 char slopet[100];
 sprintf(slopet,"%s_slope.eps",store);
 c0->Print(slopet);

 // c0->Update();

//--------------------------------------------------------------------------
// Residuals
//--------------------------------------------------------------------------
// myps.NewPage();
 
 ifstream zrese_f;
 sprintf(title,"%s/%s_Zres-calib_e.txt",txtdir,store);
 zrese_f.open(title);

 for (int i=0; i<NumBars; i++) {
    if (!zrese_f.good()) break;
    zrese_f >> dummy[0][i]
            >> zrese[0][i]
            >> dummy[1][i]
            >> zrese[1][i]
            >> dummy[2][i]
 	    >> zrese[2][i]
            >> dummy[3][i]
            >> dummy[4][i]
            >> dummy[5][i]
            >> zrese[3][i]
            >> dummy[6][i]
            >> zrese[4][i]
            >> dummy[7][i]
            >> dummy[8][i]
            >> dummy[9][i]
            >> dummy[10][i]
            >> dummy[11][i]
            >> zrese[5][i]
            >> dummy[12][i]
            >> dummy[13][i]
            >> dummy[14][i]
            >> dummy[15][i]
            >> dummy[16][i]
            >> dummy[17][i]
            >> dummy[18][i]
            >> dummy[19][i]
            >> dummy[20][i]
            >> dummy[21][i]
            >> dummy[22][i]
            >> dummy[23][i]
            >> dummy[24][i]
            >> dummy[25][i]
            >> dummy[26][i];
 }
 zrese_f.close();        
 
 ifstream zresw_f;
 sprintf(title,"%s/%s_Zres-calib_w.txt",txtdir,store);
 zresw_f.open(title);

 for (int i=0; i<NumBars; i++) {
    if (!zresw_f.good()) break;
    zresw_f >> dummy[0][i]
            >> zresw[0][i]
            >> dummy[1][i]
            >> zresw[1][i]
            >> dummy[2][i]
 	    >> zresw[2][i]
            >> dummy[3][i]
            >> dummy[4][i]
            >> dummy[5][i]
            >> zresw[3][i]
            >> dummy[6][i]
            >> zresw[4][i]
            >> dummy[7][i]
            >> dummy[8][i]
            >> dummy[9][i]
            >> dummy[10][i]
            >> dummy[11][i]
            >> zresw[5][i]
            >> dummy[12][i]
            >> dummy[13][i]
            >> dummy[14][i]
            >> dummy[15][i]
            >> dummy[16][i]
            >> dummy[17][i]
            >> dummy[18][i]
            >> dummy[19][i]
            >> dummy[20][i]
            >> dummy[21][i]
            >> dummy[22][i]
            >> dummy[23][i]
            >> dummy[24][i]
            >> dummy[25][i]
            >> dummy[26][i];
 }
 zresw_f.close();        

 for (int i=0; i<NumBars; i++) {

    if (zrese[0][i]<-2 || zrese[0][i]>2){ 
       cout << "Bar " << i << ", zrese_0 = " << zrese[0][i] << endl;
    }else{
    h1f_zres_0->Fill(zrese[0][i]); 
    }

    if (zresw[0][i]<-2 || zresw[0][i]>2){ 
       cout << "Bar " << i << ": zresw_0 = " << zresw[0][i] << endl;
    }else{
    h1f_zres_0->Fill(zresw[0][i]); 
    }

    if (zrese[1][i]<-0.006 || zrese[1][i]>0.006){ 
       cout << "Bar " << i << ", zrese_1 = " << zrese[1][i]<< endl;
    }else{
    h1f_zres_1->Fill(zrese[1][i]); 
    }

    if (zresw[1][i]<-0.006 || zresw[1][i]>0.006){ 
       cout << "Bar " << i << ": zresw_1 = " << zresw[1][i] << endl;
    }else{
    h1f_zres_1->Fill(zresw[1][i]); 
    }

    if (zrese[2][i]<-0.0009 || zrese[2][i]>0.5){ 
       cout << "Bar " << i << ", zrese_2 = " << zrese[2][i] << endl;
    }else{
    h1f_zres_2->Fill(zrese[2][i]); 
    }

    if (zresw[2][i]<-0.0009 || zresw[2][i]>0.5){ 
       cout << "Bar " << i << ": zresw_2 = " << zresw[2][i]<< endl;
    }else{
    h1f_zres_2->Fill(zresw[2][i]); 
    }

    if (zrese[3][i]<-60 || zrese[3][i]>60){ 
       cout << "Bar " << i << ", zrese_3 = " << zrese[3][i] << endl;
    }else{
    h1f_zres_3->Fill(zrese[3][i]); 
    }

    if (zresw[3][i]<-60 || zresw[3][i]>60){ 
       cout << "Bar " << i << ": zresw_3 = " << zresw[3][i] << endl;
    }else{
    h1f_zres_3->Fill(zresw[3][i]); 
    }

    if (zrese[4][i]<-0.11 || zrese[4][i]>0.09){ 
       cout << "Bar " << i << ", zrese_4 = " << zrese[4][i] << endl;
    }else{
    h1f_zres_4->Fill(zrese[4][i]); 
    }

    if (zresw[4][i]<-0.11 || zresw[4][i]>0.09 ){ 

       cout << "Bar " << i << ": zresw_4 = " << zresw[4][i] << endl;
    }else{
    h1f_zres_4->Fill(zresw[4][i]); 
    }

    if (zrese[5][i]<-500 || zrese[5][i]>500){ 
       cout << "Bar " << i << ", zrese_5 = " << zrese[5][i] << endl;
    }else{
    h1f_zres_5->Fill(zrese[5][i]); 
    }

    if (zresw[5][i]<-500 || zresw[5][i]>500){ 
       cout << "Bar " << i << ": zresw_5 = " << zresw[5][i] << endl;
    }else{
    h1f_zres_5->Fill(zresw[5][i]); 
    }
 }

 TCanvas *c5 = new TCanvas("c5","c5",800,1000);
 c5->SetFillColor(10);
 c5 ->SetBottomMargin(0.2);
 c5 ->SetLeftMargin(0.2);
 c5->Divide(2,3);

 c5->cd(1);
 Draw(h1f_zres_0,4);

 c5->cd(2);
 Draw(h1f_zres_1,4);

 c5->cd(3);
 Draw(h1f_zres_2,4);

 c5->cd(4);
 Draw(h1f_zres_3,4);

 c5->cd(5);
 Draw(h1f_zres_4,4);

 c5->cd(6);
 Draw(h1f_zres_5,4);
 char rest[100];
 sprintf(rest,"%s_res.eps",store);
 c5->Print(rest);

 // c5->Update();


//--------------------------------------------------------------------------
// Charge attenuaton
//--------------------------------------------------------------------------
// myps.NewPage();

 ifstream atte_f;
 sprintf(title,"%s/%s_Att-calib_e.txt",txtdir,store);
 atte_f.open(title);

 for (int i=0; i<NumBars; i++) {
    if (!atte_f.good()) break;
 	atte_f  >> dummy[0][i]
 		>> q_mpv_e[i]
 		>> dummy[1][i]
		>> lambda_e[i]
		>> dummy[2][i] 
		>> a_e[i]
		>> dummy[3][i] 
		>> sigma_e[i]
		>> dummy[4][i] 
        	>> dummy[5][i];
 }
 atte_f.close();      

 ifstream attw_f;
 sprintf(title,"%s/%s_Att-calib_w.txt",txtdir,store);
 attw_f.open(title);

 for (int i=0; i<NumBars; i++) {
    if (!attw_f.good()) break;
 	attw_f  >> dummy[0][i]
 		>> q_mpv_w[i]
 		>> dummy[1][i]
		>> lambda_w[i]
		>> dummy[2][i] 
		>> a_w[i]
		>> dummy[3][i] 
		>> sigma_w[i]
		>> dummy[4][i] 
        	>> dummy[5][i];
 }
 attw_f.close();      

 for (int i=0; i<NumBars; i++) {
    if (q_mpv_e[i]<0 || q_mpv_e[i]>3000){ 
       cout << "PMT " << i << " (east): q mvp e = " << q_mpv_e[i] << endl;
    }else{
       h1f_qmpv->Fill(q_mpv_e[i]); 
    }
    if (q_mpv_w[i]<0 || q_mpv_w[i]>3000){ 
       cout << "PMT " << i << " (west): q mvp w = " << q_mpv_w[i] << endl;
    }else{
       h1f_qmpv->Fill(q_mpv_w[i]); 
    }
    if (lambda_e[i]<0 || lambda_e[i]>400){ 
       cout << "PMT " << i << " (east): lambda e = " << lambda_e[i] << endl;
    }else{
       h1f_lambda->Fill(lambda_e[i]); 
    }
    if (lambda_w[i]<0 || lambda_w[i]>400){ 
       cout << "PMT " << i << " (west): lambda w = " << lambda_w[i] << endl;
    }else{
       h1f_lambda->Fill(lambda_w[i]); 
    }
    if (a_e[i]<0 || a_e[i]>0.00003){ 
       cout << "PMT " << i << " (east): a e = " << a_e[i] << endl;
    }else{
       h1f_a->Fill(a_e[i]); 
    }
    if (a_w[i]<0 || a_w[i]>0.00003){ 
       cout << "PMT " << i << " (west): a w = " << a_w[i] << endl;
    }else{
       h1f_a->Fill(a_w[i]); 
    }
    if (sigma_e[i]<0 || sigma_e[i]>500){ 
       cout << "PMT " << i << " (east): sigma e = " << sigma_e[i] << endl;
    }else{
       h1f_sigma->Fill(sigma_e[i]); 
    }
    if (sigma_w[i]<0 || sigma_w[i]>400){ 
       cout << "PMT " << i << " (west): sigma w = " << sigma_w[i] << endl;
    }else{
       h1f_sigma->Fill(sigma_w[i]); 
    }
 }

 TCanvas *c6 = new TCanvas("c6","c6",800,1000);
 c6->SetFillColor(10);
 c6 ->SetBottomMargin(0.2);
 c6 ->SetLeftMargin(0.2);
 c6->Divide(2,2);

 c6->cd(1);
 Draw(h1f_qmpv, 4);

 c6->cd(2);
 Draw(h1f_lambda, 4);

 c6->cd(3);
 Draw(h1f_a, 4);

 c6->cd(4);
 Draw(h1f_sigma, 4);
 char sigmat[100];
 sprintf(sigmat,"%s_sigma.eps",store);
 c6->Print(sigmat);

 //c6->Update();

 // myps.NewPage();

 TCanvas *c7 = new TCanvas("c7","c7",800,1000);
 c7->SetFillColor(10);
 c7 ->SetBottomMargin(0.2);
 c7 ->SetLeftMargin(0.2);
 //c7->Divide(2,2);

 for (int i=0; i<NumBars; i++) {
    h2f_resol_qmpv->Fill(q_mpv_e[i], rese[i]*1000); 
    h2f_resol_qmpv->Fill(q_mpv_w[i], resw[i]*1000); 
     //cout << "PMT " << i << " (east): resol = " << rese[i]*1000 << ", Q MPV = " 
     //                          << q_mpv_e[i] << endl;
     //cout << "PMT " << i << " (west): resol = " << resw[i]*1000 << ", Q MPV = " 
     //                         << q_mpv_w[i] << endl;
 }

 sprintf(title,"Res vs Q MPV (run range %s)", store);
 h2f_resol_qmpv->SetTitle(title);
 h2f_resol_qmpv->SetXTitle("Q MPV");
 h2f_resol_qmpv->SetYTitle("Resolution (ps)");
 h2f_resol_qmpv->GetYaxis()->SetTitleOffset(1.5);
 h2f_resol_qmpv->SetMarkerSize(3);
 h2f_resol_qmpv->Draw();

 // myps.Close();
 
 //c7->Update();
 char qtitle[100];
 sprintf(qtitle,"%s_qmip.eps", store);
 c7->Print(qtitle); 

}

//------------------------------------------------------------------------------
// Draw TH1F
//------------------------------------------------------------------------------
void Draw(TH1F *histo,        int    fColor = 0,   int fStyle = 1001,
          char *drawOpt = "", int    lColor = 1,   int lWidth = 1,
	  int   mColor  = 1,  double mSize  = 0.5, int mStyle = 20)
{
  histo->SetDirectory(0     );
  histo->SetFillColor(fColor);
  histo->SetFillStyle(fStyle);

  histo->SetLineColor  (lColor);
  histo->SetLineWidth  (lWidth);
  histo->SetMarkerColor(mColor);
  histo->SetMarkerSize (mSize );
  histo->SetMarkerStyle(mStyle);

  histo->Draw(drawOpt);
}