checkMeasured_Target.C
Go to the documentation of this file.
1 {
2 
3  bool dump = 1;
4  bool print = 1;
5 
6  TFile fxv("mX.root"); // residual in X (XZ view)
7  TFile fxz("mXz.root"); // residual in Z (XZ view)
8  TFile fyv("mY.root"); // residual in Y (YZ view)
9  TFile fyz("mYz.root"); // residual in Z (YZ view)
10 
11  TH2D *hxv = (TH2D*)((TCanvas*)fxv.Get("mX" ))->GetListOfPrimitives()->At(1);
12  TH2D *hxz = (TH2D*)((TCanvas*)fxz.Get("mXz"))->GetListOfPrimitives()->At(1);
13  TH2D *hyv = (TH2D*)((TCanvas*)fyv.Get("mY" ))->GetListOfPrimitives()->At(1);
14  TH2D *hyz = (TH2D*)((TCanvas*)fyz.Get("mYz"))->GetListOfPrimitives()->At(1);
15 
16  const int NP = 250;
17  const int NC = 100;
18  double dx[NP][NC];
19  double dy[NP][NC];
20  double dz[NP][NC];
21 
22  // open the text file holding the simulated random shifts
23  // this defines your target shifts
24  FILE *fp = fopen("/nova/app/users/gsdavies/art/align/Geometry/test/randShifts.txt","r");
25  while (!feof(fp)) {
26  int p;
27  int c;
28  double x, y, z;
29  fscanf(fp,"%d %d %lf %lf %lf\n",&p,&c,&x,&y,&z);
30  dx[p][c] = x;
31  dy[p][c] = y;
32  dz[p][c] = z;
33  }
34 
35  TH2D *goal_hxv = (TH2D*)hxv->Clone("goal_hxv");
36  TH2D *goal_hxz = (TH2D*)hxz->Clone("goal_hxz");
37  TH2D *goal_hyv = (TH2D*)hyv->Clone("goal_hyv");
38  TH2D *goal_hyz = (TH2D*)hyz->Clone("goal_hyz");
39 
40  TH2D *ratio_hxv = (TH2D*)hxv->Clone("ratio_hxv");
41  TH2D *ratio_hxz = (TH2D*)hxz->Clone("ratio_hxz");
42  TH2D *ratio_hyv = (TH2D*)hyv->Clone("ratio_hyv");
43  TH2D *ratio_hyz = (TH2D*)hyz->Clone("ratio_hyz");
44 
45  TH2D correl_xv("correl_xv",";target x-shift (x-view);measured x-shift (x-view)",100,-2.0,2.0,100,-2.0,2.0);
46  TH2D correl_xz("correl_xz",";target z-shift (x-view);measured z-shift (x-view)",100,-2.0,2.0,100,-2.0,2.0);
47  TH2D correl_yv("correl_yv",";target y-shift (y-view);measured y-shift (y-view)",100,-2.0,2.0,100,-2.0,2.0);
48  TH2D correl_yz("correl_yz",";target z-shift (y-view);measured z-shift (y-view)",100,-2.0,2.0,100,-2.0,2.0);
49 
50  correl_xv->SetStats(0);
51  correl_xz->SetStats(0);
52  correl_yv->SetStats(0);
53  correl_yz->SetStats(0);
54 
55  correl_xv->GetYaxis()->SetDecimals();
56  correl_xz->GetYaxis()->SetDecimals();
57  correl_yv->GetYaxis()->SetDecimals();
58  correl_yz->GetYaxis()->SetDecimals();
59 
60  correl_xv->GetXaxis()->SetDecimals();
61  correl_xz->GetXaxis()->SetDecimals();
62  correl_yv->GetXaxis()->SetDecimals();
63  correl_yz->GetXaxis()->SetDecimals();
64 
65  for (Int_t iP = 1; iP<=hxz->GetNbinsX(); iP++) {
66  for (Int_t iC = 1; iC<=hxz->GetNbinsY(); iC++) {
67 
68  if (hxv->GetBinContent(iP,iC)>-999) {
69 
70  //int entries = ((TH1D*)fall.Get(Form("alignment/residualmcZ_%03d_%03d",iP,iC)))->GetEntries();
71  //double rmsv = ((TH1D*)fall.Get(Form("alignment/residualmcV_%03d_%03d",iP,iC)))->GetRMS();
72  //double rmsz = ((TH1D*)fall.Get(Form("alignment/residualmcZ_%03d_%03d",iP,iC)))->GetRMS();
73 
74  goal_hxv->SetBinContent(iP,iC,-(dx[iP-1][iC-1])); // subtracting the mean here, and flipping signs
75  goal_hxz->SetBinContent(iP,iC,-(dz[iP-1][iC-1])); // subtracting the mean here, and flipping signs
76  correl_xv->Fill(-(dx[iP-1][iC-1]),hxv->GetBinContent(iP,iC));
77  correl_xz->Fill(-(dz[iP-1][iC-1]),hxz->GetBinContent(iP,iC));
78 
79  if (dump) {
80  if(hxv->GetBinContent(iP,iC)==-10.25 || hyz->GetBinContent(iP,iC)==-10.25){
81  cout << "XXX " << iP-1 << " " << iC-1 << " " <<
82  dx[iP-1][iC-1] + 0 << " " <<
83  dy[iP-1][iC-1] << " " <<
84  dz[iP-1][iC-1] + 0 << endl;
85  }
86  else{
87  cout << "XXX " << iP-1 << " " << iC-1 << " " <<
88  dx[iP-1][iC-1] + hxv->GetBinContent(iP,iC) << " " <<
89  dy[iP-1][iC-1] << " " <<
90  dz[iP-1][iC-1] + hxz->GetBinContent(iP,iC) << endl;
91  }
92  }
93  }
94  else {
95  goal_hxv->SetBinContent(iP,iC,1);
96  goal_hxz->SetBinContent(iP,iC,1);
97  }
98 
99  }
100  }
101 
102  cout << "--" << endl;
103  for (Int_t iP = 1; iP<=hyz->GetNbinsX(); iP++) {
104  for (Int_t iC = 1; iC<=hyz->GetNbinsY(); iC++) {
105  if (hyv->GetBinContent(iP,iC)>-999) {
106  goal_hyv->SetBinContent(iP,iC,-(dy[iP-1][iC-1])); // subtracting the mean here, and flipping signs
107  goal_hyz->SetBinContent(iP,iC,-(dz[iP-1][iC-1])); // subtracting the mean here, and flipping signs
108  correl_yv->Fill(-(dy[iP-1][iC-1]),hyv->GetBinContent(iP,iC));
109  correl_yz->Fill(-(dz[iP-1][iC-1]),hyz->GetBinContent(iP,iC));
110  if (dump) {
111  if(hyv->GetBinContent(iP,iC)==-10.25 || hyz->GetBinContent(iP,iC)==-10.25){
112  cout << "YYY " << iP-1 << " " << iC-1 << " " <<
113  dx[iP-1][iC-1] << " " <<
114  dy[iP-1][iC-1] + 0 << " " <<
115  dz[iP-1][iC-1] + 0 << endl;
116  }
117  else {
118  cout << "YYY " << iP-1 << " " << iC-1 << " " <<
119  dx[iP-1][iC-1] << " " <<
120  dy[iP-1][iC-1] + hyv->GetBinContent(iP,iC) << " " <<
121  dz[iP-1][iC-1] + hyz->GetBinContent(iP,iC) << endl;
122  }
123  }
124  }
125  else {
126  goal_hyv->SetBinContent(iP,iC,1);
127  goal_hyz->SetBinContent(iP,iC,1);
128  }
129  }
130  }
131 
132  ratio_hxv->Divide(hxv, goal_hxv);
133  ratio_hxz->Divide(hxz, goal_hxz);
134  ratio_hyv->Divide(hyv, goal_hyv);
135  ratio_hyz->Divide(hyz, goal_hyz);
136 
137  ratio_hxv->GetZaxis()->SetRangeUser(0,2);
138  ratio_hxz->GetZaxis()->SetRangeUser(0,2);
139  ratio_hyv->GetZaxis()->SetRangeUser(0,2);
140  ratio_hyz->GetZaxis()->SetRangeUser(0,2);
141 
142  TCanvas *canv = new TCanvas("c","c",500,500);
143 
144  if (print) {
145  correl_xv.Draw("zcol");
146  canv->Print("alignment_correl_xv_rand.png");
147  correl_xz.Draw("zcol");
148  canv->Print("alignment_correl_xz_rand.png");
149  correl_yv.Draw("zcol");
150  canv->Print("alignment_correl_yv_rand.png");
151  correl_yz.Draw("zcol");
152  canv->Print("alignment_correl_yz_rand.png");
153 
154  ratio_hxv.Draw("zcol");
155  canv->Print("alignment_ratio_hxv_rand.png");
156  ratio_hxz.Draw("zcol");
157  canv->Print("alignment_ratio_hxz_rand.png");
158  ratio_hyv.Draw("zcol");
159  canv->Print("alignment_ratio_hyv_rand.png");
160  ratio_hyz.Draw("zcol");
161  canv->Print("alignment_ratio_hyz_rand.png");
162  }
163 
164 }
TH2D * goal_hxz
TH2D * goal_hxv
TH2D correl_xv("correl_xv",";target x-shift (x-view);measured x-shift (x-view)", 100,-2.0, 2.0, 100,-2.0, 2.0)
const char * p
Definition: xmltok.h:285
TCanvas * canv
TH2D * hxv
TH2D correl_yz("correl_yz",";target z-shift (y-view);measured z-shift (y-view)", 100,-2.0, 2.0, 100,-2.0, 2.0)
TH2D correl_xz("correl_xz",";target z-shift (x-view);measured z-shift (x-view)", 100,-2.0, 2.0, 100,-2.0, 2.0)
TH2D * ratio_hyz
double dy[NP][NC]
double dx[NP][NC]
TFile fxv("mX.root")
double dz[NP][NC]
TH2D * goal_hyv
TH2D * hyz
TH2D * ratio_hxz
bool print
TFile fxz("mXz.root")
TH2D * ratio_hxv
const int NP
TFile fyz("mYz.root")
TH2D correl_yv("correl_yv",";target y-shift (y-view);measured y-shift (y-view)", 100,-2.0, 2.0, 100,-2.0, 2.0)
z
Definition: test.py:28
TH2D * ratio_hyv
OStream cout
Definition: OStream.cxx:6
const int NC
TFile fyv("mY.root")
TH2D * hxz
FILE * fp
TH2D * hyv
TH2D * goal_hyz