Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-08-05 08:15:36

0001 //step 3 with phi,r coords
0002 #include <iostream>
0003 #include <cmath>
0004 #include <vector>
0005 #include "TMath.h"
0006 #include "TVector3.h"
0007 #include "TTree.h"
0008 
0009 using namespace std;
0010 
0011 class Shifter {
0012 public:
0013 Shifter(TString sourcefilename);
0014   TFile *forward, *average;
0015   TH3F *hX, *hY, *hZ, *hR, *hPhi, *hXave, *hYave, *hZave, *hRave, *hPhiave;  
0016 };
0017 
0018 Shifter::Shifter(TString sourcefilename){
0019   //single event distortion file
0020   forward=TFile::Open(sourcefilename,"READ"); 
0021 
0022   hX=(TH3F*)forward->Get("hIntDistortionPosX");
0023   hY=(TH3F*)forward->Get("hIntDistortionPosY");
0024   hZ=(TH3F*)forward->Get("hIntDistortionPosZ");
0025 
0026   hR=(TH3F*)forward->Get("hIntDistortionPosR");
0027   hPhi=(TH3F*)forward->Get("hIntDistortionPosP");
0028 
0029   //average distortion file
0030   average=TFile::Open("/sphenix/user/rcorliss/distortion_maps/2021.04/apr07.average.real_B1.4_E-400.0.ross_phi1_sphenix_phislice_lookup_r26xp40xz40.distortion_map.hist.root","READ"); 
0031   
0032   hXave=(TH3F*)average->Get("hIntDistortionPosX");
0033   hYave=(TH3F*)average->Get("hIntDistortionPosY");
0034   hZave=(TH3F*)average->Get("hIntDistortionPosZ");
0035   
0036   hRave=(TH3F*)average->Get("hIntDistortionPosR");
0037   hPhiave=(TH3F*)average->Get("hIntDistortionPosP");
0038  
0039   //subtract average from total distortions to study fluctuations
0040   hX->Add(hXave,-1);
0041   hY->Add(hYave,-1);
0042   hZ->Add(hZave,-1);
0043   
0044   hR->Add(hRave,-1);
0045   hPhi->Add(hPhiave,-1);
0046 }
0047 
0048 int CMDistortionAnalysisPhiR(int nMaxEvents = -1) {
0049   Shifter *shifter;
0050   int nbins = 35; 
0051   double x, y, z;
0052   double low = -80.0;
0053   double high = 80.0;
0054   double deltaX, deltaY, deltaZ, deltaR, deltaPhi;
0055   int nEvents; 
0056   
0057   TCanvas *canvas=new TCanvas("canvas","CMDistortionAnalysisPhiR",2000,3000);
0058 
0059   int nsumbins = 20;
0060   int minsum = -10;
0061   int maxsum = 10;
0062   
0063   //set up summary plots
0064   TH1F *hDifferenceMeanR = new TH1F("hDifferenceMeanR", "Average Difference between R Model and True of All Events (R > 30); #Delta R (#mum)", nsumbins, minsum, maxsum);
0065     TH1F *hDifferenceStdDevR = new TH1F("hDifferenceStdDevR", "Std Dev of Difference between R Model and True of All Events (R > 30); #Delta R (#mum)", nsumbins, minsum, maxsum);
0066     
0067     TH1F *hTrueMeanR = new TH1F("hTrueMeanR", "Mean True R Distortion Model of All Events (R > 30); #Delta R (#mum)", nsumbins, minsum, maxsum);
0068     TH1F *hTrueStdDevR = new TH1F("hTrueStdDevR", "Std Dev of True R Distortion Model of All Events (R > 30); #Delta R (#mum)", nsumbins, minsum, maxsum);
0069     
0070     TH1F *hDifferenceMeanPhi = new TH1F("hDifferenceMeanPhi", "Average Difference between Phi Model and True of All Events (R > 30); #Delta Phi (#mum)", nsumbins, minsum, maxsum);
0071     TH1F *hDifferenceStdDevPhi = new TH1F("hDifferenceStdDevPhi", "Std Dev of Difference between Phi Model and True of All Events (R > 30); #Delta Phi (#mum)", nsumbins, minsum, maxsum);
0072     
0073     TH1F *hTrueMeanPhi = new TH1F("hTrueMeanPhi", "Mean True Phi Distortion Model of All Events (R > 30); #Delta Phi (#mum)", nsumbins, minsum, maxsum);
0074     TH1F *hTrueStdDevPhi = new TH1F("hTrueStdDevPhi", "Std Dev of True Phi Distortion Model of All Events (R > 30); #Delta Phi (#mum)", nsumbins, minsum, maxsum);
0075 
0076     const char * inputpattern="/sphenix/user/rcorliss/distortion_maps/2021.04/*h_Charge_*.root"; //updated
0077     
0078   //find all files that match the input string (includes wildcards)
0079   TFileCollection *filelist=new TFileCollection();
0080   filelist->Add(inputpattern);
0081   TString sourcefilename;
0082   
0083   //how many events
0084   if (nMaxEvents<0){
0085     nEvents=filelist->GetNFiles();
0086   } else if(nMaxEvents<filelist->GetNFiles()){
0087     nEvents=nMaxEvents;
0088   } else {
0089     nEvents= filelist->GetNFiles();
0090   }
0091 
0092   for (int ifile=0;ifile < nEvents;ifile++){
0093     //for each file, find all histograms in that file.
0094     sourcefilename=((TFileInfo*)(filelist->GetList()->At(ifile)))->GetCurrentUrl()->GetFile();
0095 
0096     //create shifter
0097     shifter = new Shifter(sourcefilename);
0098     
0099     TFile *plots;
0100 
0101     plots=TFile::Open(Form("CMModelsPhiR_Event%d.root",ifile),"READ");
0102 
0103     TH3F *hCartCMModelPhiR[3];
0104     hCartCMModelPhiR[0]=(TH3F*)plots->Get("hCMModelX_PhiR");
0105     hCartCMModelPhiR[1]=(TH3F*)plots->Get("hCMModelY_PhiR");
0106     hCartCMModelPhiR[2]=(TH3F*)plots->Get("hCMModelZ_PhiR");
0107 
0108     TH3F *hCylCMModelPhiR[2];
0109     hCylCMModelPhiR[0]=(TH3F*)plots->Get("hCMModelR_PhiR");
0110     hCylCMModelPhiR[1]=(TH3F*)plots->Get("hCMModelPhi_PhiR");
0111     
0112     //for forward only
0113 
0114     //same range and bins for each coordinate, binned in cm
0115     //hardcoded numbers from average distortion file's hIntDistortionPosX
0116     int nphi = 82;
0117     int nr = 54;
0118     int nz = 82;
0119     
0120     double minphi = -0.078539819;
0121     double minr = 18.884615;
0122     double minz = -1.3187500;
0123     
0124     double maxphi = 6.3617253;
0125     double maxr = 79.115387;
0126     double maxz = 106.81875;
0127 
0128     double rshiftcart, phishiftcart;
0129 
0130     int ndiff = 300;
0131     int mindiff = -20;
0132     int maxdiff = 20;
0133     
0134      TH1F *hCartesianShiftDifferencePhiR[3];
0135     hCartesianShiftDifferencePhiR[0] = new TH1F("hShiftDifferenceX_PhiR", "Difference between CM Model X and True, Phi,R binning (R > 30); #Delta X (#mum)", ndiff, mindiff, maxdiff);
0136     hCartesianShiftDifferencePhiR[1] = new TH1F("hShiftDifferenceY_PhiR", "Difference between CM Model Y and True, Phi,R binning (R > 30); #Delta Y (#mum)", ndiff, mindiff, maxdiff);
0137     hCartesianShiftDifferencePhiR[2] = new TH1F("hShiftDifferenceZ_PhiR", "Difference between CM Model Z and True, Phi,R binning (R > 30); #Delta Z (#mum)", ndiff, mindiff, maxdiff);
0138     
0139     TH1F *hCylindricalShiftDifferencePhiR[2];
0140     hCylindricalShiftDifferencePhiR[0] = new TH1F("hShiftDifferenceR_PhiR", "Difference between CM Model R and True, Phi,R binning (R > 30); #Delta R (#mum)", ndiff, mindiff, maxdiff);
0141     hCylindricalShiftDifferencePhiR[1] = new TH1F("hShiftDifferencePhi_PhiR", "Difference between CM Model Phi and True, Phi,R binning (R > 30); #Delta Phi (#mum)", ndiff, mindiff, maxdiff);
0142 
0143     TH1F *hRShiftTrue = new TH1F("hRShiftTrue", "True R Distortion Model (R > 30); #Delta R (#mum)", ndiff, mindiff, maxdiff);
0144     TH1F *hPhiShiftTrue = new TH1F("hPhiShiftTrue", "True Phi Distortion Model (R > 30); #Delta Phi (#mum)", ndiff, mindiff, maxdiff);
0145   
0146      TH2F *hCartesianDiffPhiR[6];
0147     hCartesianDiffPhiR[0] = new TH2F("hDiffXYX_PhiR", "Difference in PhiR for CM Model X, Phi,R binning; phi (rad); r (cm)",nphi,minphi,maxphi,nr,minr,maxr);
0148     hCartesianDiffPhiR[1] = new TH2F("hDiffRZX_PhiR", "Difference in RZ for CM Model X, Phi,R binning; z (cm); r (cm)", nz,minz,maxz,nr,minr,maxr);
0149     hCartesianDiffPhiR[2] = new TH2F("hDiffXYY_PhiR", "Difference in PhiR for CM Model Y, Phi,R binning; phi (rad); r (cm)",nphi,minphi,maxphi,nr,minr,maxr);
0150     hCartesianDiffPhiR[3] = new TH2F("hDiffRZY_PhiR", "Difference in RZ for CM Model Y, Phi,R binning; z (cm); r (cm)", nz,minz,maxz,nr,minr,maxr);
0151     hCartesianDiffPhiR[4] = new TH2F("hDiffXYZ_PhiR", "Difference in PhiR for CM Model Z, Phi,R binning; phi (rad); r (cm)",nphi,minphi,maxphi,nr,minr,maxr);
0152     hCartesianDiffPhiR[5] = new TH2F("hDiffRZZ_PhiR", "Difference in RZ for CM Model Z, Phi,R binning; z (cm); r (cm)", nz,minz,maxz,nr,minr,maxr);
0153     
0154     TH2F *hCylindricalDiffPhiR[4];
0155     hCylindricalDiffPhiR[0] = new TH2F("hDiffXYR_PhiR", "Difference in PhiR for CM Model R, Phi,R binning; phi (rad); r (cm)",nphi,minphi,maxphi,nr,minr,maxr);
0156     hCylindricalDiffPhiR[1] = new TH2F("hDiffRZR_PhiR", "Difference in RZ for CM Model R, Phi,R binning; z (cm); r (cm)",nz,minz,maxz,nr,minr,maxr);
0157     hCylindricalDiffPhiR[2] = new TH2F("hDiffXYPhi_PhiR", "Difference in PhiR for CM Model Phi, Phi,R binning; phi (rad); r (cm)",nphi,minphi,maxphi,nr,minr,maxr);
0158     hCylindricalDiffPhiR[3] = new TH2F("hDiffRZPhi_PhiR", "Difference in RZ for CM Model Phi, Phi,R binning; z (cm); r (cm)",nz,minz,maxz,nr,minr,maxr);
0159   
0160     TH2F *hCartesianAveDiffPhiR[6];
0161     hCartesianAveDiffPhiR[0] = new TH2F("hAveDiffXYX_PhiR", "X Model - Truth Averaged Over z, Phi,R binning (#mum); phi (rad); r (cm)",nphi,minphi,maxphi,nr,minr,maxr);
0162     hCartesianAveDiffPhiR[1] = new TH2F("hAveDiffRZX_PhiR", "X Model - Truth Averaged Over phi, Phi,R binning (#mum); z (cm); r (cm)", nz,minz,maxz,nr,minr,maxr);
0163     hCartesianAveDiffPhiR[2] = new TH2F("hAveDiffXYY_PhiR", "Y Model - Truth Averaged Over z, Phi,R binning (#mum); phi (rad); r (cm)",nphi,minphi,maxphi,nr,minr,maxr);
0164     hCartesianAveDiffPhiR[3] = new TH2F("hAveDiffRZY_PhiR", "Y Model - Truth Averaged Over phi, Phi,R binning (#mum); z (cm); r (cm)", nz,minz,maxz,nr,minr,maxr);
0165     hCartesianAveDiffPhiR[4] = new TH2F("hAveDiffXYZ_PhiR", "Z Model - Truth Averaged Over z, Phi,R binning (#mum); phi (rad); r (cm)",nphi,minphi,maxphi,nr,minr,maxr);
0166     hCartesianAveDiffPhiR[5] = new TH2F("hAveDiffRZZ_PhiR", "Z Model - Truth Averaged Over phi, Phi,R binning (#mum); z (cm); r (cm)", nz,minz,maxz,nr,minr,maxr);
0167     
0168      TH2F *hCylindricalAveDiffPhiR[4];
0169     hCylindricalAveDiffPhiR[0] = new TH2F("hAveDiffXYR_PhiR", "R Model - Truth Averaged Over z, Phi,R binning (#mum); phi (rad); r (cm)",nphi,minphi,maxphi,nr,minr,maxr);
0170     hCylindricalAveDiffPhiR[1] = new TH2F("hAveDiffRZR_PhiR", "R Model - Truth Averaged Over phi, Phi,R binning (#mum); z (cm); r (cm)", nz,minz,maxz,nr,minr,maxr);
0171     hCylindricalAveDiffPhiR[2] = new TH2F("hAveDiffXYPhi_PhiR", "Phi Model - Truth Averaged Over z, Phi,R binning (#mum); phi (rad); r (cm)",nphi,minphi,maxphi,nr,minr,maxr);
0172     hCylindricalAveDiffPhiR[3] = new TH2F("hAveDiffRZPhi_PhiR", "Phi Model - Truth Averaged Over phi, Phi,R binning (#mum); z (cm); r (cm)", nz,minz,maxz,nr,minr,maxr);
0173 
0174     TH2F *hSamplePerBinRZ = new TH2F("hSamplePerBinRZ", "Filling each rz bin; z (cm); r (cm)", nz,minz,maxz,nr,minr,maxr);
0175 
0176     TH2F *hSamplePerBinPhiR = new TH2F("hSamplePerBinPhiR", "Filling each PhiR bin; phi (rad); r (cm)",nphi,minphi,maxphi,nr,minr,maxr);
0177 
0178     TH2F *hCompareRTrue_PhiR = new TH2F("hCompareRTrue_PhiR", "Compare Difference from R Model and True, Phi,R binning (R > 30, 10 < z < 90); reco shift (#mum); true shift (#mum)",nbins,-550,550,nbins,-550,550);
0179     TH2F *hComparePhiTrue_PhiR = new TH2F("hComparePhiTrue_PhiR", "Compare Difference from Phi Model and True, Phi,R binning (R > 30, 10 < z < 90); reco shift (#mum); true shift (#mum)",nbins,-550,550,nbins,-550,550);
0180 
0181     TH2F *hRDiffvR_PhiR = new TH2F("hRDiffvR_PhiR", "Difference between R Model and True vs. r, Phi,R binning (R > 30, 10 < z < 90); r (cm); shift difference (#mum)",nr,minr,maxr,ndiff,mindiff,maxdiff);
0182     TH2F *hRDiffvZ_PhiR = new TH2F("hRDiffvZ_PhiR", "Difference between R Model and True vs. z, Phi,R binning (R > 30); z (cm); shift difference (#mum)",nz,minz,maxz,ndiff,mindiff,maxdiff);
0183     TH2F *hRDiffvPhi_PhiR = new TH2F("hRDiffvPhi_PhiR", "Difference between R Model and True vs. phi, Phi,R binning (R > 30, 10 < z < 90); phi (rad); shift difference (#mum)",nphi,minphi,maxphi,ndiff,mindiff,maxdiff);
0184 
0185     TH2F *hPhiDiffvR_PhiR = new TH2F("hPhiDiffvR_PhiR", "Difference between Phi Model and True vs. r, Phi,R binning (R > 30, 10 < z < 90); r (cm); shift difference (#mum)",nr,minr,maxr,ndiff,mindiff,maxdiff);
0186     TH2F *hPhiDiffvZ_PhiR = new TH2F("hPhiDiffvZ_PhiR", "Difference between Phi Model and True vs. z, Phi,R binning (R > 30); z (cm); shift difference (#mum)",nz,minz,maxz,ndiff,mindiff,maxdiff);
0187     TH2F *hPhiDiffvPhi_PhiR = new TH2F("hPhiDiffvPhi_PhiR", "Difference between Phi Model and True vs. phi, Phi,R binning (R > 30, 10 < z < 90); phi (rad); shift difference (#mum)",nphi,minphi,maxphi,ndiff,mindiff,maxdiff);
0188     
0189     for(int i = 1; i < nphi - 1; i++){
0190       double phi = minphi + ((maxphi - minphi)/(1.0*nphi))*(i+0.5); //center of bin
0191       for(int j = 1; j < nr - 1; j++){
0192     double r = minr + ((maxr - minr)/(1.0*nr))*(j+0.5); //center of bin
0193     for(int k = 1; k < nz - 1; k++){
0194       double z = minz + ((maxz - minz)/(1.0*nz))*(k+0.5); //center of bin
0195 
0196       double shifttrueCart[3];
0197       double shifttrueCyl[2];
0198 
0199       double shiftrecoCartPhiR[3];
0200       double differenceCartPhiR[3];
0201 
0202       double shiftrecoCylPhiR[2];
0203       double differenceCylPhiR[2];
0204 
0205       double differenceR_PhiR, differencePhi_PhiR;    
0206 
0207       int binPhiR = hCartCMModelPhiR[0]->FindBin(phi,r,z);
0208 
0209       if((r > 30.0) && (r < 76.0)){
0210         //x y and z
0211         shifttrueCart[0] = (shifter->hX->Interpolate(phi,r,z))*(1e4); //convert from cm to micron
0212         shifttrueCart[1] = (shifter->hY->Interpolate(phi,r,z))*(1e4); //convert from cm to micron 
0213         shifttrueCart[2] = (shifter->hZ->Interpolate(phi,r,z))*(1e4); //convert from cm to micron
0214         //r and phi
0215         shifttrueCyl[0] = (shifter->hR->Interpolate(phi,r,z))*(1e4); //convert from cm to micron
0216         shifttrueCyl[1] = (shifter->hPhi->Interpolate(phi,r,z))*(1e4);
0217         hRShiftTrue->Fill(shifttrueCyl[0]);
0218         hPhiShiftTrue->Fill(shifttrueCyl[1]);
0219         
0220         for(int l = 0; l < 3; l ++){
0221           shiftrecoCartPhiR[l] =  (hCartCMModelPhiR[l]->GetBinContent(binPhiR))*(1e4);
0222       
0223           differenceCartPhiR[l] = shiftrecoCartPhiR[l] - shifttrueCart[l]; 
0224 
0225           hCartesianShiftDifferencePhiR[l]->Fill(differenceCartPhiR[l]);
0226         }
0227 
0228         //r
0229         shiftrecoCylPhiR[0] =  (hCylCMModelPhiR[0]->GetBinContent(binPhiR))*(1e4);
0230         differenceCylPhiR[0] = shiftrecoCylPhiR[0] - shifttrueCyl[0]; 
0231         hCylindricalShiftDifferencePhiR[0]->Fill(differenceCylPhiR[0]);
0232           
0233         //phi
0234         shiftrecoCylPhiR[1] = r*(1e4)*(hCylCMModelPhiR[1]->GetBinContent(binPhiR));
0235         differenceCylPhiR[1] = (shiftrecoCylPhiR[1] - shifttrueCyl[1]); 
0236         hCylindricalShiftDifferencePhiR[1]->Fill(differenceCylPhiR[1]);
0237 
0238         //x
0239         hCartesianDiffPhiR[0]->Fill(phi,r, differenceCartPhiR[0]);
0240         hCartesianDiffPhiR[1]->Fill(z,r, differenceCartPhiR[0]);
0241         //y
0242         hCartesianDiffPhiR[2]->Fill(phi,r, differenceCartPhiR[1]);    
0243         hCartesianDiffPhiR[3]->Fill(z,r, differenceCartPhiR[1]);
0244         //z
0245         hCartesianDiffPhiR[4]->Fill(phi,r, differenceCartPhiR[2]);
0246         hCartesianDiffPhiR[5]->Fill(z,r, differenceCartPhiR[2]);
0247 
0248         //r
0249         hCylindricalDiffPhiR[0]->Fill(phi,r, differenceCylPhiR[0]);
0250         hCylindricalDiffPhiR[1]->Fill(z,r, differenceCylPhiR[0]);
0251 
0252         hCompareRTrue_PhiR->Fill(shiftrecoCylPhiR[0],shifttrueCyl[0]);
0253 
0254         hRDiffvR_PhiR->Fill(r,differenceCylPhiR[0],1);
0255         hRDiffvPhi_PhiR->Fill(phi,differenceCylPhiR[0],1);
0256         hRDiffvZ_PhiR->Fill(z,differenceCylPhiR[0],1);
0257 
0258         //phi 
0259         hCylindricalDiffPhiR[2]->Fill(phi,r, differenceCylPhiR[1]);
0260         hCylindricalDiffPhiR[3]->Fill(z,r, differenceCylPhiR[1]);
0261                 
0262         hComparePhiTrue_PhiR->Fill(shiftrecoCylPhiR[1],shifttrueCyl[1]);
0263         
0264         hPhiDiffvR_PhiR->Fill(r,differenceCylPhiR[1],1);
0265         hPhiDiffvPhi_PhiR->Fill(phi,differenceCylPhiR[1],1);
0266         hPhiDiffvZ_PhiR->Fill(z,differenceCylPhiR[1],1);
0267 
0268         hSamplePerBinRZ->Fill(z,r,1);
0269         hSamplePerBinPhiR->Fill(phi,r,1);
0270       }
0271     }
0272       }
0273     }
0274 
0275     //average over z
0276     for (int m = 0; m < 6; m = m+2){
0277       hCartesianAveDiffPhiR[m]->Divide(hCartesianDiffPhiR[m],hSamplePerBinPhiR);
0278     }
0279     for (int m = 0; m < 4; m = m+2){
0280       hCylindricalAveDiffPhiR[m]->Divide(hCylindricalDiffPhiR[m],hSamplePerBinPhiR);
0281     }
0282     
0283     //average over phi
0284     for (int m = 1; m < 6; m = m+2){
0285       hCartesianAveDiffPhiR[m]->Divide(hCartesianDiffPhiR[m],hSamplePerBinRZ);
0286     }
0287     for (int m = 1; m < 4; m = m+2){
0288       hCylindricalAveDiffPhiR[m]->Divide(hCylindricalDiffPhiR[m],hSamplePerBinRZ);
0289     }
0290 
0291     //summary plots
0292     hDifferenceMeanR->Fill(hCylindricalShiftDifferencePhiR[0]->GetMean(1));
0293     hDifferenceStdDevR->Fill(hCylindricalShiftDifferencePhiR[0]->GetStdDev(1));
0294 
0295     hTrueMeanR->Fill(hRShiftTrue->GetMean(1));
0296     hTrueStdDevR->Fill(hRShiftTrue->GetStdDev(1));
0297     
0298     hDifferenceMeanPhi->Fill(hCylindricalShiftDifferencePhiR[1]->GetMean(1));
0299     hDifferenceStdDevPhi->Fill(hCylindricalShiftDifferencePhiR[1]->GetStdDev(1));
0300 
0301     hTrueMeanPhi->Fill(hPhiShiftTrue->GetMean(1));
0302     hTrueStdDevPhi->Fill(hPhiShiftTrue->GetStdDev(1));
0303 
0304     for (int m = 0; m < 6; m++){
0305       hCartesianAveDiffPhiR[m]->SetStats(0);
0306     }
0307     for (int m = 0; m < 4; m++){
0308       hCylindricalAveDiffPhiR[m]->SetStats(0);
0309     }
0310   
0311 
0312     hCompareRTrue_PhiR->SetStats(0);
0313     hComparePhiTrue_PhiR->SetStats(0);
0314 
0315     hRDiffvR_PhiR->SetStats(0);
0316     hRDiffvZ_PhiR->SetStats(0);
0317     hRDiffvPhi_PhiR->SetStats(0);
0318   
0319     hPhiDiffvR_PhiR->SetStats(0);
0320     hPhiDiffvZ_PhiR->SetStats(0);
0321     hPhiDiffvPhi_PhiR->SetStats(0);
0322     
0323     TPad *c1=new TPad("c1","",0.0,0.8,1.0,0.93); //can i do an array of pads?
0324     TPad *c2=new TPad("c2","",0.0,0.64,1.0,0.77);
0325     TPad *c3=new TPad("c3","",0.0,0.48,1.0,0.61);
0326     TPad *c4=new TPad("c4","",0.0,0.32,1.0,0.45);
0327     TPad *c5=new TPad("c5","",0.0,0.16,1.0,0.29);
0328     TPad *c6=new TPad("c6","",0.0,0.0,1.0,0.13);
0329     
0330     TPad *titlepad=new TPad("titlepad","",0.0,0.96,1.0,1.0);
0331 
0332     TPad *stitlepad1=new TPad("stitlepad1","",0.0,0.93,1.0,0.96);
0333     TPad *stitlepad2=new TPad("stitlepad2","",0.0,0.77,1.0,0.8);
0334     TPad *stitlepad3=new TPad("stitlepad3","",0.0,0.61,1.0,0.64);
0335     TPad *stitlepad4=new TPad("stitlepad4","",0.0,0.45,1.0,0.48);
0336     TPad *stitlepad5=new TPad("stitlepad5","",0.0,0.29,1.0,0.32);
0337     TPad *stitlepad6=new TPad("stitlepad6","",0.0,0.13,1.0,0.16);
0338     
0339     TLatex * title = new TLatex(0.0,0.0,"");
0340 
0341     TLatex * stitle1 = new TLatex(0.0,0.0,""); //array?
0342     TLatex * stitle2 = new TLatex(0.0,0.0,"");
0343     TLatex * stitle3 = new TLatex(0.0,0.0,"");
0344     TLatex * stitle4 = new TLatex(0.0,0.0,"");
0345     TLatex * stitle5 = new TLatex(0.0,0.0,"");
0346     TLatex * stitle6 = new TLatex(0.0,0.0,"");
0347     
0348     title->SetNDC();
0349     stitle1->SetNDC();
0350     stitle2->SetNDC();
0351     stitle3->SetNDC();
0352     stitle4->SetNDC();
0353     stitle5->SetNDC();
0354     stitle6->SetNDC();
0355     
0356     title->SetTextSize(0.32);
0357     stitle1->SetTextSize(0.35);
0358     stitle2->SetTextSize(0.35);
0359     stitle3->SetTextSize(0.35);
0360     stitle4->SetTextSize(0.35);
0361     stitle5->SetTextSize(0.35);
0362     stitle6->SetTextSize(0.35);
0363     
0364     canvas->cd();
0365     c1->Draw();
0366     stitlepad1->Draw();
0367     c2->Draw();
0368     stitlepad2->Draw();
0369     c3->Draw();
0370     stitlepad3->Draw();
0371     c4->Draw();
0372     stitlepad4->Draw();
0373     c5->Draw();
0374     stitlepad5->Draw();
0375     c6->Draw();
0376     stitlepad6->Draw();
0377     titlepad->Draw();
0378 
0379     //x plots
0380     c1->Divide(4,1);
0381     c1->cd(1);
0382     hCartesianAveDiffPhiR[0]->Draw("colz");
0383     c1->cd(2);
0384     hCartesianAveDiffPhiR[1]->Draw("colz");
0385     c1->cd(3);
0386     hCartesianShiftDifferencePhiR[0]->Draw();
0387     //c1->cd(4)->Clear();  
0388     c1->cd(4);
0389     //hCMmodelSliceRvTrue->Draw("colz");
0390     hSamplePerBinRZ->Draw("colz");
0391     
0392     //y plots
0393     c2->Divide(4,1);
0394     c2->cd(1);
0395     hCartesianAveDiffPhiR[2]->Draw("colz");
0396     c2->cd(2);
0397     hCartesianAveDiffPhiR[3]->Draw("colz");
0398     c2->cd(3);
0399     hCartesianShiftDifferencePhiR[1]->Draw();
0400     //c2->cd(4)->Clear();
0401     c2->cd(4);
0402     //hStripesPerBin->Draw("colz");
0403     hSamplePerBinPhiR->Draw("colz");
0404     
0405     //r cart
0406     c3->Divide(4,1);
0407     c3->cd(1);
0408     hCylindricalAveDiffPhiR[0]->Draw("colz");
0409     c3->cd(2);
0410     hCylindricalAveDiffPhiR[1]->Draw("colz");
0411     c3->cd(3);
0412     hCylindricalShiftDifferencePhiR[0]->Draw();
0413     c3->cd(4);
0414     hRShiftTrue->Draw();
0415     
0416     //phi cart
0417     c4->Divide(4,1);
0418     c4->cd(1);
0419     hCylindricalAveDiffPhiR[2]->Draw("colz");
0420     c4->cd(2);
0421     hCylindricalAveDiffPhiR[3]->Draw("colz");
0422     c4->cd(3);
0423     hCylindricalShiftDifferencePhiR[1]->Draw();
0424     c4->cd(4);
0425     hPhiShiftTrue->Draw();
0426 
0427     //r to true comparison
0428     c5->Divide(4,1);
0429     c5->cd(1);
0430     hCompareRTrue_PhiR->Draw("colz");
0431     c5->cd(2);
0432     hRDiffvR_PhiR->Draw("colz");
0433     c5->cd(3);
0434     hRDiffvZ_PhiR->Draw("colz");
0435     c5->cd(4);
0436     hRDiffvPhi_PhiR->Draw("colz");
0437 
0438     //phi to true comparison
0439     c6->Divide(4,1);
0440     c6->cd(1);
0441     hComparePhiTrue_PhiR->Draw("colz");
0442     c6->cd(2);
0443     hPhiDiffvR_PhiR->Draw("colz");
0444     c6->cd(3);
0445     hPhiDiffvZ_PhiR->Draw("colz");
0446     c6->cd(4);
0447     hPhiDiffvPhi_PhiR->Draw("colz");
0448 
0449     titlepad->cd();
0450     titlepad->Clear();
0451     title->DrawLatex(0.01,0.4,Form("Event %d; %s", ifile, sourcefilename.Data())); 
0452     title->Draw();
0453     
0454     stitlepad1->cd();
0455     stitlepad1->Clear();
0456     stitle1->DrawLatex(0.45,0.2,"X Model"); 
0457     stitle1->Draw();
0458      
0459     stitlepad2->cd();
0460     stitlepad2->Clear();
0461     stitle2->DrawLatex(0.45,0.2,"Y Model"); 
0462     stitle2->Draw();
0463 
0464     stitlepad3->cd();
0465     stitlepad3->Clear();
0466     stitle3->DrawLatex(0.45,0.2,"R Model"); 
0467     stitle3->Draw();
0468 
0469     stitlepad4->cd();
0470     stitlepad4->Clear();
0471     stitle4->DrawLatex(0.45,0.2,"Phi Model"); 
0472     stitle4->Draw();
0473 
0474     stitlepad5->cd();
0475     stitlepad5->Clear();
0476     stitle5->DrawLatex(0.4,0.2,"Comparing R Model to True"); 
0477     stitle5->Draw();
0478 
0479     stitlepad6->cd();
0480     stitlepad6->Clear();
0481     stitle6->DrawLatex(0.4,0.2,"Comparing Phi Model to True"); 
0482     stitle6->Draw();
0483 
0484     if(ifile == 0){ 
0485       //if(ifile == 1){
0486       canvas->Print("CMDistortionAnalysisPhiR.pdf(","pdf");
0487     } else if((ifile == 1) || (ifile == nEvents - 1)){
0488       canvas->Print("CMDistortionAnalysisPhiR.pdf","pdf");
0489     }
0490   }
0491 
0492   TCanvas *summary = new TCanvas("summary","ShiftPlotsSummary",2000,3000);
0493 
0494   TPad *sumtitlepad = new TPad("sumtitlepad","",0.0,0.96,1.0,1.0);
0495   TPad *sumplots = new TPad("sumplotspad","",0.0,0.0,1.0,0.96);
0496 
0497   TLatex *sumtitle = new TLatex(0.0,0.0,"");
0498 
0499   sumtitle->SetNDC();
0500   sumtitle->SetTextSize(0.4);
0501 
0502   summary->cd();
0503   sumplots->Draw();
0504   sumtitlepad->Draw();
0505 
0506   sumplots->Divide(4,6);
0507   sumplots->cd(1);
0508   hDifferenceMeanR->Draw();
0509   sumplots->cd(2);
0510   hDifferenceStdDevR->Draw();
0511   sumplots->cd(3);
0512   hTrueMeanR->Draw();
0513   sumplots->cd(4);
0514   hTrueStdDevR->Draw();
0515   sumplots->cd(5);
0516   hDifferenceMeanPhi->Draw();
0517   sumplots->cd(6);
0518   hDifferenceStdDevPhi->Draw();
0519   sumplots->cd(7);
0520   hTrueMeanPhi->Draw();
0521   sumplots->cd(8);
0522   hTrueStdDevPhi->Draw();
0523   sumplots->cd(9);
0524   sumplots->cd(10)->Clear();
0525   sumplots->cd(11)->Clear();
0526   sumplots->cd(12)->Clear();
0527   sumplots->cd(13)->Clear();
0528   sumplots->cd(14)->Clear();
0529   sumplots->cd(15)->Clear();
0530   sumplots->cd(16)->Clear();
0531   sumplots->cd(17)->Clear();
0532   sumplots->cd(18)->Clear();
0533   sumplots->cd(19)->Clear();
0534   sumplots->cd(20)->Clear();
0535   sumplots->cd(21)->Clear();
0536   sumplots->cd(22)->Clear();
0537   sumplots->cd(23)->Clear();
0538   sumplots->cd(24)->Clear();
0539 
0540   sumtitlepad->cd();
0541   sumtitlepad->Clear();
0542   sumtitle->DrawLatex(0.4,0.4,"Summary of Events"); 
0543   summary->Print("CMDistortionAnalysisPhiR.pdf)","pdf");
0544 
0545   return 0;
0546 }