ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
CMDistortionAnalysisPhiR.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file CMDistortionAnalysisPhiR.C
1 //step 3 with phi,r coords
2 #include <iostream>
3 #include <cmath>
4 #include <vector>
5 #include "TMath.h"
6 #include "TVector3.h"
7 #include "TTree.h"
8 
9 using namespace std;
10 
11 class Shifter {
12 public:
13 Shifter(TString sourcefilename);
14  TFile *forward, *average;
15  TH3F *hX, *hY, *hZ, *hR, *hPhi, *hXave, *hYave, *hZave, *hRave, *hPhiave;
16 };
17 
18 Shifter::Shifter(TString sourcefilename){
19  //single event distortion file
20  forward=TFile::Open(sourcefilename,"READ");
21 
22  hX=(TH3F*)forward->Get("hIntDistortionPosX");
23  hY=(TH3F*)forward->Get("hIntDistortionPosY");
24  hZ=(TH3F*)forward->Get("hIntDistortionPosZ");
25 
26  hR=(TH3F*)forward->Get("hIntDistortionPosR");
27  hPhi=(TH3F*)forward->Get("hIntDistortionPosP");
28 
29  //average distortion file
30  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");
31 
32  hXave=(TH3F*)average->Get("hIntDistortionPosX");
33  hYave=(TH3F*)average->Get("hIntDistortionPosY");
34  hZave=(TH3F*)average->Get("hIntDistortionPosZ");
35 
36  hRave=(TH3F*)average->Get("hIntDistortionPosR");
37  hPhiave=(TH3F*)average->Get("hIntDistortionPosP");
38 
39  //subtract average from total distortions to study fluctuations
40  hX->Add(hXave,-1);
41  hY->Add(hYave,-1);
42  hZ->Add(hZave,-1);
43 
44  hR->Add(hRave,-1);
45  hPhi->Add(hPhiave,-1);
46 }
47 
48 int CMDistortionAnalysisPhiR(int nMaxEvents = -1) {
49  Shifter *shifter;
50  int nbins = 35;
51  double x, y, z;
52  double low = -80.0;
53  double high = 80.0;
54  double deltaX, deltaY, deltaZ, deltaR, deltaPhi;
55  int nEvents;
56 
57  TCanvas *canvas=new TCanvas("canvas","CMDistortionAnalysisPhiR",2000,3000);
58 
59  int nsumbins = 20;
60  int minsum = -10;
61  int maxsum = 10;
62 
63  //set up summary plots
64  TH1F *hDifferenceMeanR = new TH1F("hDifferenceMeanR", "Average Difference between R Model and True of All Events (R > 30); #Delta R (#mum)", nsumbins, minsum, maxsum);
65  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);
66 
67  TH1F *hTrueMeanR = new TH1F("hTrueMeanR", "Mean True R Distortion Model of All Events (R > 30); #Delta R (#mum)", nsumbins, minsum, maxsum);
68  TH1F *hTrueStdDevR = new TH1F("hTrueStdDevR", "Std Dev of True R Distortion Model of All Events (R > 30); #Delta R (#mum)", nsumbins, minsum, maxsum);
69 
70  TH1F *hDifferenceMeanPhi = new TH1F("hDifferenceMeanPhi", "Average Difference between Phi Model and True of All Events (R > 30); #Delta Phi (#mum)", nsumbins, minsum, maxsum);
71  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);
72 
73  TH1F *hTrueMeanPhi = new TH1F("hTrueMeanPhi", "Mean True Phi Distortion Model of All Events (R > 30); #Delta Phi (#mum)", nsumbins, minsum, maxsum);
74  TH1F *hTrueStdDevPhi = new TH1F("hTrueStdDevPhi", "Std Dev of True Phi Distortion Model of All Events (R > 30); #Delta Phi (#mum)", nsumbins, minsum, maxsum);
75 
76  const char * inputpattern="/sphenix/user/rcorliss/distortion_maps/2021.04/*h_Charge_*.root"; //updated
77 
78  //find all files that match the input string (includes wildcards)
79  TFileCollection *filelist=new TFileCollection();
80  filelist->Add(inputpattern);
81  TString sourcefilename;
82 
83  //how many events
84  if (nMaxEvents<0){
85  nEvents=filelist->GetNFiles();
86  } else if(nMaxEvents<filelist->GetNFiles()){
87  nEvents=nMaxEvents;
88  } else {
89  nEvents= filelist->GetNFiles();
90  }
91 
92  for (int ifile=0;ifile < nEvents;ifile++){
93  //for each file, find all histograms in that file.
94  sourcefilename=((TFileInfo*)(filelist->GetList()->At(ifile)))->GetCurrentUrl()->GetFile();
95 
96  //create shifter
97  shifter = new Shifter(sourcefilename);
98 
99  TFile *plots;
100 
101  plots=TFile::Open(Form("CMModelsPhiR_Event%d.root",ifile),"READ");
102 
103  TH3F *hCartCMModelPhiR[3];
104  hCartCMModelPhiR[0]=(TH3F*)plots->Get("hCMModelX_PhiR");
105  hCartCMModelPhiR[1]=(TH3F*)plots->Get("hCMModelY_PhiR");
106  hCartCMModelPhiR[2]=(TH3F*)plots->Get("hCMModelZ_PhiR");
107 
108  TH3F *hCylCMModelPhiR[2];
109  hCylCMModelPhiR[0]=(TH3F*)plots->Get("hCMModelR_PhiR");
110  hCylCMModelPhiR[1]=(TH3F*)plots->Get("hCMModelPhi_PhiR");
111 
112  //for forward only
113 
114  //same range and bins for each coordinate, binned in cm
115  //hardcoded numbers from average distortion file's hIntDistortionPosX
116  int nphi = 82;
117  int nr = 54;
118  int nz = 82;
119 
120  double minphi = -0.078539819;
121  double minr = 18.884615;
122  double minz = -1.3187500;
123 
124  double maxphi = 6.3617253;
125  double maxr = 79.115387;
126  double maxz = 106.81875;
127 
128  double rshiftcart, phishiftcart;
129 
130  int ndiff = 300;
131  int mindiff = -20;
132  int maxdiff = 20;
133 
134  TH1F *hCartesianShiftDifferencePhiR[3];
135  hCartesianShiftDifferencePhiR[0] = new TH1F("hShiftDifferenceX_PhiR", "Difference between CM Model X and True, Phi,R binning (R > 30); #Delta X (#mum)", ndiff, mindiff, maxdiff);
136  hCartesianShiftDifferencePhiR[1] = new TH1F("hShiftDifferenceY_PhiR", "Difference between CM Model Y and True, Phi,R binning (R > 30); #Delta Y (#mum)", ndiff, mindiff, maxdiff);
137  hCartesianShiftDifferencePhiR[2] = new TH1F("hShiftDifferenceZ_PhiR", "Difference between CM Model Z and True, Phi,R binning (R > 30); #Delta Z (#mum)", ndiff, mindiff, maxdiff);
138 
139  TH1F *hCylindricalShiftDifferencePhiR[2];
140  hCylindricalShiftDifferencePhiR[0] = new TH1F("hShiftDifferenceR_PhiR", "Difference between CM Model R and True, Phi,R binning (R > 30); #Delta R (#mum)", ndiff, mindiff, maxdiff);
141  hCylindricalShiftDifferencePhiR[1] = new TH1F("hShiftDifferencePhi_PhiR", "Difference between CM Model Phi and True, Phi,R binning (R > 30); #Delta Phi (#mum)", ndiff, mindiff, maxdiff);
142 
143  TH1F *hRShiftTrue = new TH1F("hRShiftTrue", "True R Distortion Model (R > 30); #Delta R (#mum)", ndiff, mindiff, maxdiff);
144  TH1F *hPhiShiftTrue = new TH1F("hPhiShiftTrue", "True Phi Distortion Model (R > 30); #Delta Phi (#mum)", ndiff, mindiff, maxdiff);
145 
146  TH2F *hCartesianDiffPhiR[6];
147  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);
148  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);
149  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);
150  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);
151  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);
152  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);
153 
154  TH2F *hCylindricalDiffPhiR[4];
155  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);
156  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);
157  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);
158  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);
159 
160  TH2F *hCartesianAveDiffPhiR[6];
161  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);
162  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);
163  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);
164  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);
165  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);
166  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);
167 
168  TH2F *hCylindricalAveDiffPhiR[4];
169  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);
170  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);
171  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);
172  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);
173 
174  TH2F *hSamplePerBinRZ = new TH2F("hSamplePerBinRZ", "Filling each rz bin; z (cm); r (cm)", nz,minz,maxz,nr,minr,maxr);
175 
176  TH2F *hSamplePerBinPhiR = new TH2F("hSamplePerBinPhiR", "Filling each PhiR bin; phi (rad); r (cm)",nphi,minphi,maxphi,nr,minr,maxr);
177 
178  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);
179  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);
180 
181  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);
182  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);
183  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);
184 
185  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);
186  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);
187  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);
188 
189  for(int i = 1; i < nphi - 1; i++){
190  double phi = minphi + ((maxphi - minphi)/(1.0*nphi))*(i+0.5); //center of bin
191  for(int j = 1; j < nr - 1; j++){
192  double r = minr + ((maxr - minr)/(1.0*nr))*(j+0.5); //center of bin
193  for(int k = 1; k < nz - 1; k++){
194  double z = minz + ((maxz - minz)/(1.0*nz))*(k+0.5); //center of bin
195 
196  double shifttrueCart[3];
197  double shifttrueCyl[2];
198 
199  double shiftrecoCartPhiR[3];
200  double differenceCartPhiR[3];
201 
202  double shiftrecoCylPhiR[2];
203  double differenceCylPhiR[2];
204 
205  double differenceR_PhiR, differencePhi_PhiR;
206 
207  int binPhiR = hCartCMModelPhiR[0]->FindBin(phi,r,z);
208 
209  if((r > 30.0) && (r < 76.0)){
210  //x y and z
211  shifttrueCart[0] = (shifter->hX->Interpolate(phi,r,z))*(1e4); //convert from cm to micron
212  shifttrueCart[1] = (shifter->hY->Interpolate(phi,r,z))*(1e4); //convert from cm to micron
213  shifttrueCart[2] = (shifter->hZ->Interpolate(phi,r,z))*(1e4); //convert from cm to micron
214  //r and phi
215  shifttrueCyl[0] = (shifter->hR->Interpolate(phi,r,z))*(1e4); //convert from cm to micron
216  shifttrueCyl[1] = (shifter->hPhi->Interpolate(phi,r,z))*(1e4);
217  hRShiftTrue->Fill(shifttrueCyl[0]);
218  hPhiShiftTrue->Fill(shifttrueCyl[1]);
219 
220  for(int l = 0; l < 3; l ++){
221  shiftrecoCartPhiR[l] = (hCartCMModelPhiR[l]->GetBinContent(binPhiR))*(1e4);
222 
223  differenceCartPhiR[l] = shiftrecoCartPhiR[l] - shifttrueCart[l];
224 
225  hCartesianShiftDifferencePhiR[l]->Fill(differenceCartPhiR[l]);
226  }
227 
228  //r
229  shiftrecoCylPhiR[0] = (hCylCMModelPhiR[0]->GetBinContent(binPhiR))*(1e4);
230  differenceCylPhiR[0] = shiftrecoCylPhiR[0] - shifttrueCyl[0];
231  hCylindricalShiftDifferencePhiR[0]->Fill(differenceCylPhiR[0]);
232 
233  //phi
234  shiftrecoCylPhiR[1] = r*(1e4)*(hCylCMModelPhiR[1]->GetBinContent(binPhiR));
235  differenceCylPhiR[1] = (shiftrecoCylPhiR[1] - shifttrueCyl[1]);
236  hCylindricalShiftDifferencePhiR[1]->Fill(differenceCylPhiR[1]);
237 
238  //x
239  hCartesianDiffPhiR[0]->Fill(phi,r, differenceCartPhiR[0]);
240  hCartesianDiffPhiR[1]->Fill(z,r, differenceCartPhiR[0]);
241  //y
242  hCartesianDiffPhiR[2]->Fill(phi,r, differenceCartPhiR[1]);
243  hCartesianDiffPhiR[3]->Fill(z,r, differenceCartPhiR[1]);
244  //z
245  hCartesianDiffPhiR[4]->Fill(phi,r, differenceCartPhiR[2]);
246  hCartesianDiffPhiR[5]->Fill(z,r, differenceCartPhiR[2]);
247 
248  //r
249  hCylindricalDiffPhiR[0]->Fill(phi,r, differenceCylPhiR[0]);
250  hCylindricalDiffPhiR[1]->Fill(z,r, differenceCylPhiR[0]);
251 
252  hCompareRTrue_PhiR->Fill(shiftrecoCylPhiR[0],shifttrueCyl[0]);
253 
254  hRDiffvR_PhiR->Fill(r,differenceCylPhiR[0],1);
255  hRDiffvPhi_PhiR->Fill(phi,differenceCylPhiR[0],1);
256  hRDiffvZ_PhiR->Fill(z,differenceCylPhiR[0],1);
257 
258  //phi
259  hCylindricalDiffPhiR[2]->Fill(phi,r, differenceCylPhiR[1]);
260  hCylindricalDiffPhiR[3]->Fill(z,r, differenceCylPhiR[1]);
261 
262  hComparePhiTrue_PhiR->Fill(shiftrecoCylPhiR[1],shifttrueCyl[1]);
263 
264  hPhiDiffvR_PhiR->Fill(r,differenceCylPhiR[1],1);
265  hPhiDiffvPhi_PhiR->Fill(phi,differenceCylPhiR[1],1);
266  hPhiDiffvZ_PhiR->Fill(z,differenceCylPhiR[1],1);
267 
268  hSamplePerBinRZ->Fill(z,r,1);
269  hSamplePerBinPhiR->Fill(phi,r,1);
270  }
271  }
272  }
273  }
274 
275  //average over z
276  for (int m = 0; m < 6; m = m+2){
277  hCartesianAveDiffPhiR[m]->Divide(hCartesianDiffPhiR[m],hSamplePerBinPhiR);
278  }
279  for (int m = 0; m < 4; m = m+2){
280  hCylindricalAveDiffPhiR[m]->Divide(hCylindricalDiffPhiR[m],hSamplePerBinPhiR);
281  }
282 
283  //average over phi
284  for (int m = 1; m < 6; m = m+2){
285  hCartesianAveDiffPhiR[m]->Divide(hCartesianDiffPhiR[m],hSamplePerBinRZ);
286  }
287  for (int m = 1; m < 4; m = m+2){
288  hCylindricalAveDiffPhiR[m]->Divide(hCylindricalDiffPhiR[m],hSamplePerBinRZ);
289  }
290 
291  //summary plots
292  hDifferenceMeanR->Fill(hCylindricalShiftDifferencePhiR[0]->GetMean(1));
293  hDifferenceStdDevR->Fill(hCylindricalShiftDifferencePhiR[0]->GetStdDev(1));
294 
295  hTrueMeanR->Fill(hRShiftTrue->GetMean(1));
296  hTrueStdDevR->Fill(hRShiftTrue->GetStdDev(1));
297 
298  hDifferenceMeanPhi->Fill(hCylindricalShiftDifferencePhiR[1]->GetMean(1));
299  hDifferenceStdDevPhi->Fill(hCylindricalShiftDifferencePhiR[1]->GetStdDev(1));
300 
301  hTrueMeanPhi->Fill(hPhiShiftTrue->GetMean(1));
302  hTrueStdDevPhi->Fill(hPhiShiftTrue->GetStdDev(1));
303 
304  for (int m = 0; m < 6; m++){
305  hCartesianAveDiffPhiR[m]->SetStats(0);
306  }
307  for (int m = 0; m < 4; m++){
308  hCylindricalAveDiffPhiR[m]->SetStats(0);
309  }
310 
311 
312  hCompareRTrue_PhiR->SetStats(0);
313  hComparePhiTrue_PhiR->SetStats(0);
314 
315  hRDiffvR_PhiR->SetStats(0);
316  hRDiffvZ_PhiR->SetStats(0);
317  hRDiffvPhi_PhiR->SetStats(0);
318 
319  hPhiDiffvR_PhiR->SetStats(0);
320  hPhiDiffvZ_PhiR->SetStats(0);
321  hPhiDiffvPhi_PhiR->SetStats(0);
322 
323  TPad *c1=new TPad("c1","",0.0,0.8,1.0,0.93); //can i do an array of pads?
324  TPad *c2=new TPad("c2","",0.0,0.64,1.0,0.77);
325  TPad *c3=new TPad("c3","",0.0,0.48,1.0,0.61);
326  TPad *c4=new TPad("c4","",0.0,0.32,1.0,0.45);
327  TPad *c5=new TPad("c5","",0.0,0.16,1.0,0.29);
328  TPad *c6=new TPad("c6","",0.0,0.0,1.0,0.13);
329 
330  TPad *titlepad=new TPad("titlepad","",0.0,0.96,1.0,1.0);
331 
332  TPad *stitlepad1=new TPad("stitlepad1","",0.0,0.93,1.0,0.96);
333  TPad *stitlepad2=new TPad("stitlepad2","",0.0,0.77,1.0,0.8);
334  TPad *stitlepad3=new TPad("stitlepad3","",0.0,0.61,1.0,0.64);
335  TPad *stitlepad4=new TPad("stitlepad4","",0.0,0.45,1.0,0.48);
336  TPad *stitlepad5=new TPad("stitlepad5","",0.0,0.29,1.0,0.32);
337  TPad *stitlepad6=new TPad("stitlepad6","",0.0,0.13,1.0,0.16);
338 
339  TLatex * title = new TLatex(0.0,0.0,"");
340 
341  TLatex * stitle1 = new TLatex(0.0,0.0,""); //array?
342  TLatex * stitle2 = new TLatex(0.0,0.0,"");
343  TLatex * stitle3 = new TLatex(0.0,0.0,"");
344  TLatex * stitle4 = new TLatex(0.0,0.0,"");
345  TLatex * stitle5 = new TLatex(0.0,0.0,"");
346  TLatex * stitle6 = new TLatex(0.0,0.0,"");
347 
348  title->SetNDC();
349  stitle1->SetNDC();
350  stitle2->SetNDC();
351  stitle3->SetNDC();
352  stitle4->SetNDC();
353  stitle5->SetNDC();
354  stitle6->SetNDC();
355 
356  title->SetTextSize(0.32);
357  stitle1->SetTextSize(0.35);
358  stitle2->SetTextSize(0.35);
359  stitle3->SetTextSize(0.35);
360  stitle4->SetTextSize(0.35);
361  stitle5->SetTextSize(0.35);
362  stitle6->SetTextSize(0.35);
363 
364  canvas->cd();
365  c1->Draw();
366  stitlepad1->Draw();
367  c2->Draw();
368  stitlepad2->Draw();
369  c3->Draw();
370  stitlepad3->Draw();
371  c4->Draw();
372  stitlepad4->Draw();
373  c5->Draw();
374  stitlepad5->Draw();
375  c6->Draw();
376  stitlepad6->Draw();
377  titlepad->Draw();
378 
379  //x plots
380  c1->Divide(4,1);
381  c1->cd(1);
382  hCartesianAveDiffPhiR[0]->Draw("colz");
383  c1->cd(2);
384  hCartesianAveDiffPhiR[1]->Draw("colz");
385  c1->cd(3);
386  hCartesianShiftDifferencePhiR[0]->Draw();
387  //c1->cd(4)->Clear();
388  c1->cd(4);
389  //hCMmodelSliceRvTrue->Draw("colz");
390  hSamplePerBinRZ->Draw("colz");
391 
392  //y plots
393  c2->Divide(4,1);
394  c2->cd(1);
395  hCartesianAveDiffPhiR[2]->Draw("colz");
396  c2->cd(2);
397  hCartesianAveDiffPhiR[3]->Draw("colz");
398  c2->cd(3);
399  hCartesianShiftDifferencePhiR[1]->Draw();
400  //c2->cd(4)->Clear();
401  c2->cd(4);
402  //hStripesPerBin->Draw("colz");
403  hSamplePerBinPhiR->Draw("colz");
404 
405  //r cart
406  c3->Divide(4,1);
407  c3->cd(1);
408  hCylindricalAveDiffPhiR[0]->Draw("colz");
409  c3->cd(2);
410  hCylindricalAveDiffPhiR[1]->Draw("colz");
411  c3->cd(3);
412  hCylindricalShiftDifferencePhiR[0]->Draw();
413  c3->cd(4);
414  hRShiftTrue->Draw();
415 
416  //phi cart
417  c4->Divide(4,1);
418  c4->cd(1);
419  hCylindricalAveDiffPhiR[2]->Draw("colz");
420  c4->cd(2);
421  hCylindricalAveDiffPhiR[3]->Draw("colz");
422  c4->cd(3);
423  hCylindricalShiftDifferencePhiR[1]->Draw();
424  c4->cd(4);
425  hPhiShiftTrue->Draw();
426 
427  //r to true comparison
428  c5->Divide(4,1);
429  c5->cd(1);
430  hCompareRTrue_PhiR->Draw("colz");
431  c5->cd(2);
432  hRDiffvR_PhiR->Draw("colz");
433  c5->cd(3);
434  hRDiffvZ_PhiR->Draw("colz");
435  c5->cd(4);
436  hRDiffvPhi_PhiR->Draw("colz");
437 
438  //phi to true comparison
439  c6->Divide(4,1);
440  c6->cd(1);
441  hComparePhiTrue_PhiR->Draw("colz");
442  c6->cd(2);
443  hPhiDiffvR_PhiR->Draw("colz");
444  c6->cd(3);
445  hPhiDiffvZ_PhiR->Draw("colz");
446  c6->cd(4);
447  hPhiDiffvPhi_PhiR->Draw("colz");
448 
449  titlepad->cd();
450  titlepad->Clear();
451  title->DrawLatex(0.01,0.4,Form("Event %d; %s", ifile, sourcefilename.Data()));
452  title->Draw();
453 
454  stitlepad1->cd();
455  stitlepad1->Clear();
456  stitle1->DrawLatex(0.45,0.2,"X Model");
457  stitle1->Draw();
458 
459  stitlepad2->cd();
460  stitlepad2->Clear();
461  stitle2->DrawLatex(0.45,0.2,"Y Model");
462  stitle2->Draw();
463 
464  stitlepad3->cd();
465  stitlepad3->Clear();
466  stitle3->DrawLatex(0.45,0.2,"R Model");
467  stitle3->Draw();
468 
469  stitlepad4->cd();
470  stitlepad4->Clear();
471  stitle4->DrawLatex(0.45,0.2,"Phi Model");
472  stitle4->Draw();
473 
474  stitlepad5->cd();
475  stitlepad5->Clear();
476  stitle5->DrawLatex(0.4,0.2,"Comparing R Model to True");
477  stitle5->Draw();
478 
479  stitlepad6->cd();
480  stitlepad6->Clear();
481  stitle6->DrawLatex(0.4,0.2,"Comparing Phi Model to True");
482  stitle6->Draw();
483 
484  if(ifile == 0){
485  //if(ifile == 1){
486  canvas->Print("CMDistortionAnalysisPhiR.pdf(","pdf");
487  } else if((ifile == 1) || (ifile == nEvents - 1)){
488  canvas->Print("CMDistortionAnalysisPhiR.pdf","pdf");
489  }
490  }
491 
492  TCanvas *summary = new TCanvas("summary","ShiftPlotsSummary",2000,3000);
493 
494  TPad *sumtitlepad = new TPad("sumtitlepad","",0.0,0.96,1.0,1.0);
495  TPad *sumplots = new TPad("sumplotspad","",0.0,0.0,1.0,0.96);
496 
497  TLatex *sumtitle = new TLatex(0.0,0.0,"");
498 
499  sumtitle->SetNDC();
500  sumtitle->SetTextSize(0.4);
501 
502  summary->cd();
503  sumplots->Draw();
504  sumtitlepad->Draw();
505 
506  sumplots->Divide(4,6);
507  sumplots->cd(1);
508  hDifferenceMeanR->Draw();
509  sumplots->cd(2);
510  hDifferenceStdDevR->Draw();
511  sumplots->cd(3);
512  hTrueMeanR->Draw();
513  sumplots->cd(4);
514  hTrueStdDevR->Draw();
515  sumplots->cd(5);
516  hDifferenceMeanPhi->Draw();
517  sumplots->cd(6);
518  hDifferenceStdDevPhi->Draw();
519  sumplots->cd(7);
520  hTrueMeanPhi->Draw();
521  sumplots->cd(8);
522  hTrueStdDevPhi->Draw();
523  sumplots->cd(9);
524  sumplots->cd(10)->Clear();
525  sumplots->cd(11)->Clear();
526  sumplots->cd(12)->Clear();
527  sumplots->cd(13)->Clear();
528  sumplots->cd(14)->Clear();
529  sumplots->cd(15)->Clear();
530  sumplots->cd(16)->Clear();
531  sumplots->cd(17)->Clear();
532  sumplots->cd(18)->Clear();
533  sumplots->cd(19)->Clear();
534  sumplots->cd(20)->Clear();
535  sumplots->cd(21)->Clear();
536  sumplots->cd(22)->Clear();
537  sumplots->cd(23)->Clear();
538  sumplots->cd(24)->Clear();
539 
540  sumtitlepad->cd();
541  sumtitlepad->Clear();
542  sumtitle->DrawLatex(0.4,0.4,"Summary of Events");
543  summary->Print("CMDistortionAnalysisPhiR.pdf)","pdf");
544 
545  return 0;
546 }