Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-08-05 08:14:43

0001 static const int N_ENERGIES = 8;
0002 static const int N_CUTLEVELS = 8;
0003 #include "load_files.C"
0004 void plot_resolution() {
0005 
0006     gStyle->SetOptFit(1111);
0007     gStyle->SetOptTitle(1);
0008     gStyle->SetStatW(0.3);
0009     
0010     int scan = 2;
0011     int location = 1;
0012     int energies[N_ENERGIES];
0013     
0014     string filenames[N_ENERGIES];
0015     
0016     load_files(scan,location,energies,filenames);
0017 
0018     TFile *fin[N_ENERGIES];
0019     TTree *trees[N_ENERGIES];
0020 
0021     float Etotal_t, Emax_t, E3by3_t, E5by5_t, C1_t,  C2_inner_t,  C2_outer_t, C2_inner_new_t, C2_outer_new_t;
0022     float Horz_HODO_R0_t, Horz_HODO_R1_t, Horz_HODO_R2_t, Horz_HODO_R3_t, Horz_HODO_R4_t;
0023     float Horz_HODO_R5_t, Horz_HODO_R6_t, Horz_HODO_R7_t, Vert_HODO_R0_t, Vert_HODO_R1_t;
0024     float Vert_HODO_R2_t, Vert_HODO_R3_t, Vert_HODO_R4_t, Vert_HODO_R5_t, Vert_HODO_R6_t, Vert_HODO_R7_t;
0025     float Veto1_t, Veto2_t, Veto3_t, Veto4_t, TowerID_t,Tower_column_t,Tower_row_t;
0026     float TowerE_column_0_t;
0027     float TowerE_column_1_t;
0028     float TowerE_column_2_t;
0029     float TowerE_column_3_t;
0030     float TowerE_column_4_t;
0031     float TowerE_column_5_t;
0032     float TowerE_column_6_t;
0033     float TowerE_column_7_t;
0034     float TowerE_row_0_t;
0035     float TowerE_row_1_t;
0036     float TowerE_row_2_t;
0037     float TowerE_row_3_t;
0038     float TowerE_row_4_t;
0039     float TowerE_row_5_t;
0040     float TowerE_row_6_t;
0041     float TowerE_row_7_t;
0042     float HorzTowerID_t;
0043     float VertTowerID_t;
0044     float SaveHoriz_TowerID0_t;
0045     float SaveHoriz_TowerID1_t;
0046     float SaveHoriz_TowerID2_t;
0047     float SaveHoriz_TowerID3_t;
0048     float SaveHoriz_TowerID4_t;
0049     float SaveHoriz_TowerID5_t;
0050     float SaveHoriz_TowerID6_t;
0051     float SaveHoriz_TowerID7_t;
0052 
0053     float SaveVert_TowerID0_t;
0054     float SaveVert_TowerID1_t;
0055     float SaveVert_TowerID2_t;
0056     float SaveVert_TowerID3_t;
0057     float SaveVert_TowerID4_t;
0058     float SaveVert_TowerID5_t;
0059     float SaveVert_TowerID6_t;
0060     float SaveVert_TowerID7_t;
0061 
0062     for(int i=0; i<N_ENERGIES; i++){
0063          fin[i] = new TFile(filenames[i].c_str());
0064 
0065          ostringstream name;
0066          name << "test1";
0067          trees[i] = (TTree *)fin[i]->Get(name.str().c_str());
0068 
0069          name.str("");
0070          name.clear();
0071 
0072          name << "tree_" << energies[i];
0073          trees[i]->SetName(name.str().c_str());
0074          trees[i]->SetBranchAddress("Etotal_t",&Etotal_t);
0075          trees[i]->SetBranchAddress("E3by3_t",&E3by3_t);
0076          trees[i]->SetBranchAddress("Emax_t",&Emax_t);
0077 
0078          trees[i]->SetBranchAddress("E5by5_t",&E5by5_t);
0079          trees[i]->SetBranchAddress("C1_t",&C1_t);
0080          trees[i]->SetBranchAddress("C2_inner_t",&C2_inner_t);
0081          trees[i]->SetBranchAddress("C2_outer_t",&C2_outer_t);
0082 
0083         if (scan==3) {
0084             trees[i]->SetBranchAddress("C2_inner_new_t",&C2_inner_new_t);
0085             trees[i]->SetBranchAddress("C2_outer_new_t",&C2_outer_new_t);
0086         }
0087  
0088          trees[i]->SetBranchAddress("Horz_HODO_R0_t",&Horz_HODO_R0_t);
0089          trees[i]->SetBranchAddress("Horz_HODO_R1_t",&Horz_HODO_R1_t);
0090          trees[i]->SetBranchAddress("Horz_HODO_R2_t",&Horz_HODO_R2_t);
0091          trees[i]->SetBranchAddress("Horz_HODO_R3_t",&Horz_HODO_R3_t);
0092          trees[i]->SetBranchAddress("Horz_HODO_R4_t",&Horz_HODO_R4_t);
0093          trees[i]->SetBranchAddress("Horz_HODO_R5_t",&Horz_HODO_R5_t);
0094          trees[i]->SetBranchAddress("Horz_HODO_R6_t",&Horz_HODO_R6_t);
0095          trees[i]->SetBranchAddress("Horz_HODO_R7_t",&Horz_HODO_R7_t);
0096          trees[i]->SetBranchAddress("Vert_HODO_R0_t",&Vert_HODO_R0_t);
0097          trees[i]->SetBranchAddress("Vert_HODO_R1_t",&Vert_HODO_R1_t);
0098          trees[i]->SetBranchAddress("Vert_HODO_R2_t",&Vert_HODO_R2_t);
0099          trees[i]->SetBranchAddress("Vert_HODO_R3_t",&Vert_HODO_R3_t);
0100          trees[i]->SetBranchAddress("Vert_HODO_R4_t",&Vert_HODO_R4_t);
0101          trees[i]->SetBranchAddress("Vert_HODO_R5_t",&Vert_HODO_R5_t);
0102          trees[i]->SetBranchAddress("Vert_HODO_R6_t",&Vert_HODO_R6_t);
0103          trees[i]->SetBranchAddress("Vert_HODO_R7_t",&Vert_HODO_R7_t);
0104          trees[i]->SetBranchAddress("Veto1_t",&Veto1_t);
0105          trees[i]->SetBranchAddress("Veto2_t",&Veto2_t);
0106          trees[i]->SetBranchAddress("Veto3_t",&Veto3_t);
0107          trees[i]->SetBranchAddress("Veto4_t",&Veto4_t);
0108          trees[i]->SetBranchAddress("TowerID_t",&TowerID_t);
0109          trees[i]->SetBranchAddress("Tower_column_t",&Tower_column_t);
0110          trees[i]->SetBranchAddress("Tower_row_t",&Tower_row_t);
0111          trees[i]->SetBranchAddress("TowerE_column_0_t",&TowerE_column_0_t);
0112          trees[i]->SetBranchAddress("TowerE_column_1_t",&TowerE_column_1_t);
0113          trees[i]->SetBranchAddress("TowerE_column_2_t",&TowerE_column_2_t);
0114          trees[i]->SetBranchAddress("TowerE_column_3_t",&TowerE_column_3_t);
0115          trees[i]->SetBranchAddress("TowerE_column_4_t",&TowerE_column_4_t);
0116          trees[i]->SetBranchAddress("TowerE_column_5_t",&TowerE_column_5_t);
0117          trees[i]->SetBranchAddress("TowerE_column_6_t",&TowerE_column_6_t);
0118          trees[i]->SetBranchAddress("TowerE_column_7_t",&TowerE_column_7_t);
0119          trees[i]->SetBranchAddress("TowerE_row_0_t",&TowerE_row_0_t);
0120          trees[i]->SetBranchAddress("TowerE_row_1_t",&TowerE_row_1_t);
0121          trees[i]->SetBranchAddress("TowerE_row_2_t",&TowerE_row_2_t);
0122          trees[i]->SetBranchAddress("TowerE_row_3_t",&TowerE_row_3_t);
0123          trees[i]->SetBranchAddress("TowerE_row_4_t",&TowerE_row_4_t);
0124          trees[i]->SetBranchAddress("TowerE_row_5_t",&TowerE_row_5_t);
0125          trees[i]->SetBranchAddress("TowerE_row_6_t",&TowerE_row_6_t);
0126          trees[i]->SetBranchAddress("TowerE_row_7_t",&TowerE_row_7_t);
0127          trees[i]->SetBranchAddress("HorzTowerID_t",&HorzTowerID_t);
0128          trees[i]->SetBranchAddress("VertTowerID_t",&VertTowerID_t);
0129 
0130          trees[i]->SetBranchAddress("SaveHoriz_TowerID0_t",&SaveHoriz_TowerID0_t);
0131          trees[i]->SetBranchAddress("SaveHoriz_TowerID1_t",&SaveHoriz_TowerID1_t);
0132          trees[i]->SetBranchAddress("SaveHoriz_TowerID2_t",&SaveHoriz_TowerID2_t);
0133          trees[i]->SetBranchAddress("SaveHoriz_TowerID3_t",&SaveHoriz_TowerID3_t);
0134          trees[i]->SetBranchAddress("SaveHoriz_TowerID4_t",&SaveHoriz_TowerID4_t);
0135          trees[i]->SetBranchAddress("SaveHoriz_TowerID5_t",&SaveHoriz_TowerID5_t);
0136          trees[i]->SetBranchAddress("SaveHoriz_TowerID6_t",&SaveHoriz_TowerID6_t);
0137          trees[i]->SetBranchAddress("SaveHoriz_TowerID7_t",&SaveHoriz_TowerID7_t);
0138          
0139          trees[i]->SetBranchAddress("SaveVert_TowerID0_t",&SaveVert_TowerID0_t);
0140          trees[i]->SetBranchAddress("SaveVert_TowerID1_t",&SaveVert_TowerID1_t);
0141          trees[i]->SetBranchAddress("SaveVert_TowerID2_t",&SaveVert_TowerID2_t);
0142          trees[i]->SetBranchAddress("SaveVert_TowerID3_t",&SaveVert_TowerID3_t);
0143          trees[i]->SetBranchAddress("SaveVert_TowerID4_t",&SaveVert_TowerID4_t);
0144          trees[i]->SetBranchAddress("SaveVert_TowerID5_t",&SaveVert_TowerID5_t);
0145          trees[i]->SetBranchAddress("SaveVert_TowerID6_t",&SaveVert_TowerID6_t);
0146          trees[i]->SetBranchAddress("SaveVert_TowerID7_t",&SaveVert_TowerID7_t);
0147      }
0148 
0149     TF1 *fgauss[N_CUTLEVELS][N_ENERGIES];  //for the fit arrays
0150     TH1D *hE5x5[N_CUTLEVELS][N_ENERGIES]; //for the histogram arrays
0151      
0152     double range_histogram[2][N_ENERGIES];
0153     range_histogram[0][0] = 0;
0154     range_histogram[1][0] = 5;
0155     range_histogram[0][1] = 0;
0156     range_histogram[1][1] = 15;
0157     range_histogram[0][2] = 0;
0158     range_histogram[1][2] = 15;
0159     range_histogram[0][3] = 0;
0160     range_histogram[1][3] = 18;
0161     range_histogram[0][4] = 1;
0162     range_histogram[1][4] = 22;
0163     range_histogram[0][5] = 1;
0164     range_histogram[1][5] = 30;
0165     range_histogram[0][6] = 1;
0166     range_histogram[1][6] = 22;
0167     range_histogram[0][7] = 1;
0168     range_histogram[1][7] = 25;
0169 
0170     double range_fit[2][N_ENERGIES];
0171     //scan 3
0172     if (scan==3) {
0173         range_fit[0][0] = 1.7;
0174         range_fit[1][0] = 2.7;
0175         range_fit[0][1] = 3.5;
0176         range_fit[1][1] = 5.5;
0177         range_fit[0][2] = 5.7;
0178         range_fit[1][2] = 7.8;
0179         range_fit[0][3] = 8.0;
0180         range_fit[1][3] = 10.8;
0181         range_fit[0][4] = 11.5;
0182         range_fit[1][4] = 15.8;
0183         range_fit[0][5] = 15.0;
0184         range_fit[1][5] = 20;
0185         range_fit[0][6] = 10.2;
0186         range_fit[1][6] = 13.4;
0187         range_fit[0][7] = 11.6;
0188         range_fit[1][7] = 16.0;
0189     }
0190     
0191     if (scan==2) {
0192         if (location==1) {
0193             range_fit[0][0] = 0.7;
0194             range_fit[1][0] = 1.4;
0195             range_fit[0][1] = 1.5;
0196             range_fit[1][1] = 2.7;
0197             range_fit[0][2] = 2.5;
0198             range_fit[1][2] = 4.0;
0199             range_fit[0][3] = 3.3;
0200             range_fit[1][3] = 4.9;
0201             range_fit[0][4] = 5.2;
0202             range_fit[1][4] = 7.5;
0203             range_fit[0][5] = 7.3;
0204             range_fit[1][5] = 9.9;
0205             range_fit[0][6] = 10.9;
0206             range_fit[1][6] = 15.0;
0207             range_fit[0][7] = 14.7;
0208             range_fit[1][7] = 19.2;
0209 
0210         }
0211         if (location==3) {
0212             range_fit[0][0] = 1.8;
0213             range_fit[1][0] = 3.5;
0214             range_fit[0][1] = 1.8;
0215             range_fit[1][1] = 2.8;
0216             range_fit[0][2] = 2.8;
0217             range_fit[1][2] = 4.5;
0218             range_fit[0][3] = 4.0;
0219             range_fit[1][3] = 5.8;
0220             range_fit[0][4] = 6.5;
0221             range_fit[1][4] = 8.0;
0222             range_fit[0][5] = 8.0;
0223             range_fit[1][5] = 11.5;
0224             range_fit[0][6] = 13.5;
0225             range_fit[1][6] = 16.5;
0226             range_fit[0][7] = 17.0;
0227             range_fit[1][7] = 21.5;
0228         }
0229     }
0230     
0231      /*
0232         cutlevel 0: no cuts
0233         cutlevel 1: cherenkov
0234         cutlevel 2: veto
0235         cutlevel 3: hodo
0236         cutlevel 4: cherenkov + veto
0237         cutlevel 5: cherenkov + hodo
0238         cutlevel 6: veto + hodo
0239         cutlevel 7: cherenkov + veto + hodo
0240       */
0241     string cutlevel[N_CUTLEVELS];
0242     cutlevel[0] = "no_cuts";
0243     cutlevel[1] = "cherenkov";
0244     cutlevel[2] = "veto";
0245     cutlevel[3] = "hodoscope";
0246     cutlevel[4] = "cherenkov+veto";
0247     cutlevel[5] = "cherenkov+hodo";
0248     cutlevel[6] = "veto+hodo";
0249     cutlevel[7] = "cherenkov+veto+hodo";
0250 
0251     for(int i=0; i<N_ENERGIES; i++){
0252 
0253         cout << "starting: " << energies[i] << " GeV" << endl;
0254 
0255         //setup the energy distribution plots
0256         int entries = trees[i]->GetEntries();
0257         for(int j=0; j<N_CUTLEVELS; j++){
0258             //name the plots
0259             ostringstream name;
0260             name << "hE5x5_" << cutlevel[j] << "_" << energies[i]<< " GeV";
0261             hE5x5[j][i] = new TH1D(name.str().c_str(),name.str().c_str(),200,range_histogram[0][i],range_histogram[1][i]);
0262             hE5x5[j][i]->Sumw2();
0263 
0264             name.str("");
0265             name.clear();
0266             name << "fgauss_" << j << "_" << i;
0267             fgauss[j][i] = new TF1(name.str().c_str(),"gaus",range_fit[0][i],range_fit[1][i]);
0268         }
0269     
0270         for(j = 0; j<entries; j++){//loop to fill histograms
0271             trees[i]->GetEntry(j);
0272 
0273             bool passCherenkov = false;
0274             bool passVeto = false;
0275             //no cuts
0276             hE5x5[0][i]->Fill(E5by5_t);
0277           
0278             //Cherenkov cut
0279             if (scan==3) {
0280                 if(fabs(C2_inner_new_t + C2_outer_new_t) > 100){
0281                     passCherenkov = true;
0282                     hE5x5[1][i]->Fill(E5by5_t);
0283                 }
0284             }
0285             if (scan==2) {
0286                 if(fabs(C2_inner_t + C2_outer_t) > 100){
0287                     passCherenkov = true;
0288                     hE5x5[1][i]->Fill(E5by5_t);
0289                 }
0290             }
0291             //veto cut
0292             if (Veto1_t < 15 || Veto2_t < 15 || Veto3_t < 15 || Veto4_t < 15) {
0293                 passVeto =  true;
0294                 hE5x5[2][i]->Fill(E5by5_t);
0295             }
0296 
0297             //vertical hodoscope cut
0298             bool passVertical = false;
0299             if (scan==2) {
0300                 if (location==1) {
0301                     if (fabs(Vert_HODO_R2_t) > 30||fabs(Vert_HODO_R3_t)>30) {
0302                         passVertical = true;
0303                     }
0304                 }
0305                 if (location==2) {
0306                     if (fabs(Vert_HODO_R2_t) > 30||fabs(Vert_HODO_R3_t)>30) {
0307                         passVertical = true;
0308                     }
0309                 }
0310                 if (location==3) {
0311                     if (fabs(Vert_HODO_R4_t) > 30) {
0312                         passVertical = true;
0313                     }
0314                 }
0315             }
0316             if (scan==3) {
0317                 if (fabs(Vert_HODO_R5_t) > 30||fabs(Vert_HODO_R6_t)>30) {
0318                     passVertical = true;
0319                 }
0320             }
0321             //horizontal hodoscope cut
0322             bool passHorizontal = false;
0323             if (scan==2) {
0324                 if (location==1) {
0325                     if (fabs(Horz_HODO_R5_t) > 30) {
0326                         passHorizontal = true;
0327                     }
0328                 }
0329                 if (location==2) {
0330                     if (fabs(Horz_HODO_R3_t) > 30)||fabs(Horz_HODO_R4_t) > 30) {
0331                         passHorizontal = true;
0332                     }
0333                 }
0334                 if (location==3) {
0335                     if (fabs(Horz_HODO_R3_t) > 30)||fabs(Horz_HODO_R4_t) > 30) {
0336                         passHorizontal = true;
0337                     }
0338                 }
0339             }
0340             if (scan==3) {
0341                 if (fabs(Horz_HODO_R4_t) > 30) {
0342                     passHorizontal = true;
0343                 }
0344             }
0345             
0346             if (passHorizontal && passVertical) {
0347                 hE5x5[3][i]->Fill(E5by5_t);
0348             }
0349 
0350             //veto+chrenkov cut
0351             if(passCherenkov && passVeto){
0352                 hE5x5[4][i]->Fill(E5by5_t);
0353             }
0354             //cherenkov+hodo
0355             if (passVertical && passHorizontal && passCherenkov) {
0356                 hE5x5[5][i]->Fill(E5by5_t);
0357             }
0358             //veto+hodo
0359             if (passVeto && passVertical && passHorizontal) {
0360                 hE5x5[6][i]->Fill(E5by5_t);
0361             }
0362             //veto + hodo +cherenkov
0363             if(passCherenkov && passVeto && passVertical && passHorizontal) {
0364                 hE5x5[7][i]->Fill(E5by5_t);
0365             }
0366         }
0367     }
0368 
0369     TCanvas *c2[N_CUTLEVELS];  //canvas with energy resolutions
0370     
0371     //~~~
0372     for(int i=0; i<N_CUTLEVELS; i++){
0373         ostringstream name;
0374         if (i==0)
0375             name << "no cuts_";
0376         if (i==1)
0377             name<< "cherenkov cuts ";
0378         if (i==2)
0379             name<< "veto cuts ";
0380         if (i==3)
0381             name<<" hodoscope cuts";
0382         if (i==4)
0383             name<< "chrenkov+veto";
0384         if (i==5)
0385             name<<"cherenkov+hodoscopes";
0386         if (i==6)
0387             name<<"veto+hodoscopes";
0388         if (i==7)
0389             name<<"cherenkov + veto + hodoscopes";
0390 
0391         c2[i] = new TCanvas(name.str().c_str(),name.str().c_str(),1500,700);//for the histograms of energies, all the cuts
0392         c2[i]->Divide(4,2);
0393         //define the naming convaention on the overall canvas
0394         for(int j =0; j<N_ENERGIES; j++)//loop through the enrgies on the same plot
0395         {
0396             c2[i]->cd(j+1);
0397             gPad->SetLogy();
0398             hE5x5[i][j]->GetXaxis()->SetTitle("energy in 5x5 (GeV) ");
0399             hE5x5[i][j]->Draw();
0400             hE5x5[i][j]->Fit(fgauss[i][j],"R");
0401 
0402         }//end j loop through 8 energy plots
0403         
0404     }//end of i loop
0405 
0406      //==========================================================================================
0407     TCanvas *c3 = new TCanvas("c3","EMCAL Energy Resolution",1500,700);
0408     c3->Divide(4,2);
0409     
0410     Float_t *mean [N_ENERGIES];
0411     
0412     for (int j=0; j<N_CUTLEVELS; j++)
0413     {
0414         Float_t  getMEAN[N_ENERGIES];
0415         Float_t getSIGMA[N_ENERGIES];
0416         Float_t getSIGMAerror[N_ENERGIES];
0417         Float_t getMEANerror[N_ENERGIES];
0418         
0419         Float_t getAllerror[N_ENERGIES];
0420         Float_t y_1[N_ENERGIES];  //stores the ratio of sigma/mean
0421         Float_t y_2[N_ENERGIES];  //stores the mean
0422         
0423         cout<<cutlevel[j]<<":"<<endl;
0424         for (int i=0;i<N_ENERGIES;i++)
0425         {
0426             getMEAN[i]=fgauss[j][i]->GetParameter(1);
0427             getSIGMA[i]=fgauss[j][i]->GetParameter(2);
0428             getMEANerror[i]=fgauss[j][i]->GetParError(1);
0429             getSIGMAerror[i]=fgauss[j][i]->GetParError(2);
0430             
0431             getAllerror[i]=(getSIGMA[i]/getMEAN[i])*TMath::Sqrt( (getSIGMAerror[i]/getSIGMA[i])*(getSIGMAerror[i]/getSIGMA[i])+(getMEANerror[i]/getMEAN[i])*(getMEANerror[i]/getMEAN[i]) );
0432             y_1[i]=getSIGMA[i]/getMEAN[i];
0433             y_2[i]=getMEAN[i];
0434             cout<<"All error "<<i<<" is :"<<getAllerror[i]<<endl;
0435             cout<<"All y "<<i<<" is :"<<y_1[i]<<endl;
0436         
0437         }  //end of energy loop
0438         c3->cd(j+1);
0439         c3->SetFillColor(0);
0440         c3->SetGrid();
0441         c3->GetFrame()->SetBorderSize(12);
0442         
0443         const Int_t n1=8;
0444         TF1 *f2=new TF1("f2","TMath::Sqrt([0]*[0]+([1]*[1])*TMath::Power(x,-1))",0.5,28);  //convolved
0445         if (scan==3) {
0446             Float_t x1[n1]  = {2,4,6,8,12,16,24,28};
0447         }
0448         if (scan==2) {
0449             Float_t x1[n1]  = {1,2,3,4,6,8,12,16};
0450         }
0451         Float_t ex1[n1] = {0,0,0,0,0,0,0,0};
0452         TGraphErrors *grEtotal = new TGraphErrors(n1,x1,y_1,ex1,getAllerror);
0453         grEtotal->Fit("f2","R");//fit the Etotal convolved
0454         grEtotal->SetMarkerColor(kBlack);
0455         grEtotal->SetMarkerStyle(20);
0456         ostringstream name;
0457         name<<cutlevel[j];
0458         grEtotal->SetTitle(name.str().c_str());
0459         grEtotal->GetXaxis()->SetTitle("Energy (GeV)");
0460         grEtotal->GetXaxis()->SetRangeUser(0,20);
0461         grEtotal->GetYaxis()->SetTitle("#sigma(E)/<E>");
0462         grEtotal->GetHistogram()->SetMaximum(0.25);    // along
0463         grEtotal->GetHistogram()->SetMinimum(0.02);   //   Y
0464         grEtotal->Draw("AP");
0465     }//end of N_cutlevels
0466     //c3->SaveAs("EMCal_Energy_Resolution_scan3.pdf");
0467 
0468     TCanvas *c6 = new TCanvas("c6","EMCAL Linearity & Energy Resolution ",1500,700);
0469     
0470     c6->Divide(2,1);
0471     c6->cd(2);
0472     const Int_t n1=8;
0473     TF1 *f2=new TF1("f2","TMath::Sqrt([0]*[0]+([1]*[1])*TMath::Power(x,-1))",0.5,28);  //convolved
0474     if (scan==3) {
0475         Float_t x1[n1]  = {2,4,6,8,12,16,24,28};
0476     }
0477     if (scan==2) {
0478         Float_t x1[n1]  = {1,2,3,4,6,8,12,16};
0479     }
0480     Float_t ex1[n1] = {0,0,0,0,0,0,0,0};
0481     TGraphErrors *grEtotal = new TGraphErrors(n1,x1,y_1,ex1,getAllerror);
0482     grEtotal->Fit("f2","R");//fit the Etotal convolved
0483     grEtotal->SetMarkerColor(kBlack);
0484     grEtotal->SetMarkerStyle(20);
0485     ostringstream name;
0486     name<<cutlevel[j];
0487     grEtotal->SetTitle("cherenkov+veto+hodo");
0488     grEtotal->GetXaxis()->SetTitle("Energy (GeV)");
0489     grEtotal->GetXaxis()->SetRangeUser(0,30);
0490     grEtotal->GetYaxis()->SetTitle("#sigma(E)/<E>");
0491     grEtotal->GetHistogram()->SetMaximum(0.3);    // along
0492     grEtotal->GetHistogram()->SetMinimum(-0.1);   //   Y
0493     grEtotal->Draw("AP");
0494     c6->cd(1);
0495     //make a linearity plot
0496     gStyle->SetStatX(0.9);
0497     gStyle->SetStatY(0.35);
0498 
0499     TF1 * f_calo_l = new TF1("f_calo_l", "pol1", 0.5, 20);
0500     TGraphErrors *grLinearity = new TGraphErrors(n1,x1,y_2,ex1,getAllerror);
0501     grLinearity->Fit("f_calo_l","R");//fit the Etotal convolved
0502     grLinearity->SetMarkerColor(kBlack);
0503     grLinearity->SetMarkerStyle(20);
0504     grLinearity->SetTitle("Linearity 3rd scan");
0505     grLinearity->GetXaxis()->SetTitle("Input Energy (GeV)");
0506     grLinearity->GetXaxis()->SetRangeUser(1,30);
0507     grLinearity->GetYaxis()->SetTitle("Measured Energy (GeV)");
0508     grLinearity->Draw("AP");
0509     TGraphErrors *y_x = new TGraphErrors(n1,x1,x1,ex1,ex1);
0510     y_x->SetLineColor(kGray);
0511     y_x->Draw("same");
0512     
0513 }
0514 
0515