StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
alignTPC.C
1 void alignTPC(){
2 
3  gSystem->Load("libPhysics");
4  gStyle->SetFrameBorderMode(0);
5  gStyle->SetPadBorderMode(0);
6 
7  TVector3 EAO[25], EBO[25], EDO[25];
8  TVector3 WAO[25], WBO[25], WDO[25];
9 
10  TVector3 tempVect(0,0,0);
11  TVector3 crossVect(0,0,0);
12  TVector3 ortVect(0,0,0);
13  TVector3 uniVect(0,0,0);
14  TVector3 nullVect(0,0,0);
15 
16  TVector3 axisSurveyX(1.,0.,0.);
17  TVector3 axisSurveyY(0.,1.,0.);
18  TVector3 axisSurveyZ(0.,0.,1.);
19 
20  double tAngle =0.;
21 
22  EAO[13].SetXYZ(-0.27465, -2.31703, 0.47066);
23  EAO[14].SetXYZ(-0.47431, -2.31705, 0.27091);
24  EAO[15].SetXYZ(-0.54726, -2.31713, -0.00182);
25  EAO[16].SetXYZ(-0.47394, -2.31718, -0.27460);
26  EAO[17].SetXYZ(-0.27426, -2.31711, -0.47420);
27  EAO[18].SetXYZ(-0.00134, -2.31699, -0.54718);
28  EAO[19].SetXYZ(0.27141, -2.31675, -0.47399);
29  EAO[20].SetXYZ(0.47104, -2.31663, -0.27427);
30  EAO[21].SetXYZ(0.54403, -2.31668, -0.00137);
31  EAO[22].SetXYZ(0.47082, -2.31674, 0.27141);
32  EAO[23].SetXYZ(0.27101, -2.31683, 0.47091);
33  EAO[24].SetXYZ(-0.00182, -2.31695, 0.54387);
34 
35  EBO[13].SetXYZ( -0.59588, -2.31725, 1.02658);
36  EBO[14].SetXYZ(-1.03053, -2.31719, 0.59175);
37  EBO[15].SetXYZ(-1.18926, -2.31731, -0.00211);
38  EBO[16].SetXYZ(-1.03017, -2.31692, -0.59588);
39  EBO[17].SetXYZ(-0.59496, -2.31720, -1.03037);
40  EBO[18].SetXYZ(-0.00110, -2.31688, -1.18911);
41  EBO[19].SetXYZ(0., 0., 0.);
42  EBO[20].SetXYZ(1.02724, -2.31631, -0.59508);
43  EBO[21].SetXYZ(1.18602, -2.31632, -0.00115);
44  EBO[22].SetXYZ(1.02673, -2.31650, 0.59253);
45  EBO[23].SetXYZ(0.59184, -2.31677, 1.02697);
46  EBO[24].SetXYZ(-0.00208, -2.31696, 1.18566);
47 
48  EDO[13].SetXYZ(-0.95525, -2.31701, 1.64841);
49  EDO[14].SetXYZ(0., 0., 0.);
50  EDO[15].SetXYZ(0., 0., 0.);
51  EDO[16].SetXYZ(-1.65192, -2.31734, -0.95513);
52  EDO[17].SetXYZ(-0.95400, -2.31674, -1.65246);
53  EDO[18].SetXYZ(0., 0., 0.);
54  EDO[19].SetXYZ(0.95194, -2.31608, -1.65172);
55  EDO[20].SetXYZ(1.64927, -2.31585, -0.95392);
56  EDO[21].SetXYZ(0., 0., 0.);
57  EDO[22].SetXYZ(0., 0., 0.);
58  EDO[23].SetXYZ(0.95070, -2.31622, 1.64905);
59  EDO[24].SetXYZ(-0.00234, -2.31674, 1.90405);
60 
61  WAO[1].SetXYZ(0., 0., 0.);
62  WAO[2].SetXYZ(0.46923, 2.31308, 0.27104);
63  WAO[3].SetXYZ(0.54225, 2.31310, -0.00181);
64  WAO[4].SetXYZ(0.46903, 2.31310, -0.27453);
65  WAO[5].SetXYZ(0.26930, 2.31289, -0.47419);
66  WAO[6].SetXYZ(-0.00342, 2.31277, -0.54732);
67  WAO[7].SetXYZ(-0.27622, 2.31275, -0.47430);
68  WAO[8].SetXYZ(-0.47585, 2.31267, -0.27457);
69  WAO[9].SetXYZ(-0.54889, 2.31278, -0.00186);
70  WAO[10].SetXYZ(-0.47583, 2.31282, 0.27090);
71  WAO[11].SetXYZ(-0.27613, 2.31298, 0.47050);
72  WAO[12].SetXYZ(-0.00335, 2.31298, 0.54363);
73 
74  WBO[1].SetXYZ(0.59054, 2.31298, 1.02659);
75  WBO[2].SetXYZ(1.02529, 2.31314, 0.59191);
76  WBO[3].SetXYZ(1.18429, 2.31314, -0.00186);
77  WBO[4].SetXYZ(1.02510, 2.31299, -0.59558);
78  WBO[5].SetXYZ(0.59031, 2.31274, -1.03022);
79  WBO[6].SetXYZ(-0.00346, 2.31251, -1.18937);
80  WBO[7].SetXYZ(-0.59727, 2.31251, -1.03034);
81  WBO[8].SetXYZ(-1.03189, 2.31226, -0.59552);
82  WBO[9].SetXYZ(-1.19092, 2.31227, -0.00178);
83  WBO[10].SetXYZ(-1.03191, 2.31243, 0.59190);
84  WBO[11].SetXYZ(-0.59707, 2.31270, 1.02648);
85  WBO[12].SetXYZ(-0.00330, 2.31279, 1.18557);
86 
87  WDO[1].SetXYZ(0.94962, 2.31289, 1.64844);
88  WDO[2].SetXYZ(1.64729, 2.31318, 0.95101);
89  WDO[3].SetXYZ(1.90240, 2.31331, -0.00176);
90  WDO[4].SetXYZ(0., 0., 0.);
91  WDO[5].SetXYZ(0.94945, 2.31280, -1.65202);
92  WDO[6].SetXYZ(-0.00344, 2.31223, -1.90736);
93  WDO[7].SetXYZ(-0.95633, 2.31209, -1.65220);
94  WDO[8].SetXYZ(0., 0., 0.);
95  WDO[9].SetXYZ(0., 0., 0.);
96  WDO[10].SetXYZ(-1.65380, 2.31169, 0.95098);
97  WDO[11].SetXYZ(0., 0., 0.);
98  WDO[12].SetXYZ(0., 0., 0.);
99 
100  TH3F *tpcCenterInSurvey = new TH3F("tpcCenterInSurvey","tpcCenterInSurvey",101,-10.1,10.1,101,-10.1,10.1,101,-10.1,10.1);
101  TH1F *tpcCenterInSurveyX = new TH1F("tpcCenterInSurveyX","tpcCenterInSurveyX",201,-10.05,10.05);
102  tpcCenterInSurveyX->SetXTitle("X_{survey} (mm)");
103  TH1F *tpcCenterInSurveyY = new TH1F("tpcCenterInSurveyY","tpcCenterInSurveyY",201,-10.05,10.05);
104  tpcCenterInSurveyY->SetXTitle("Y_{survey} (mm)");
105  TH1F *tpcCenterInSurveyZ = new TH1F("tpcCenterInSurveyZ","tpcCenterInSurveyZ",201,-10.05,10.05);
106  tpcCenterInSurveyZ->SetXTitle("Z_{survey} (mm)");
107 
108  TH3F *tpcWTiltInSurvey = new TH3F("tpcWTiltInSurvey","tpcWTiltInSurvey",101,-0.00101,0.00101,101,-0.00101,0.00101,101,-0.00101,0.00101);
109  TH1F *tpcWTiltInSurveyX = new TH1F("tpcWTiltInSurveyX","tpcWTiltInSurveyX",201,-0.001005,0.001005);
110  TH1F *tpcWTiltInSurveyY = new TH1F("tpcWTiltInSurveyY","tpcWTiltInSurveyY",201,-0.001005,0.001005);
111  TH1F *tpcWTiltInSurveyZ = new TH1F("tpcWTiltInSurveyZ","tpcWTiltInSurveyZ",201,-0.001005,0.001005);
112 
113  TH3F *tpcETiltInSurvey = new TH3F("tpcETiltInSurvey","tpcETiltInSurvey",101,-1010,1010,101,-1010,1010,101,-1010,1010);
114  TH1F *tpcETiltInSurveyX = new TH1F("tpcETiltInSurveyX","tpcETiltInSurveyX",201,-1005,1005);
115  TH1F *tpcETiltInSurveyY = new TH1F("tpcETiltInSurveyY","tpcETiltInSurveyY",201,-1005,1005);
116  TH1F *tpcETiltInSurveyZ = new TH1F("tpcETiltInSurveyZ","tpcETiltInSurveyZ",201,-1005,1005);
117 
118  TH1F *tpcWXInSurveyXY = new TH1F("tpcWXInSurveyXY","tpcWXInSurveyXY",201,-1.005,1.005);
119  tpcWXInSurveyXY->SetXTitle("Angle (mrad)");
120  TH1F *tpcWXInSurveyXZ = new TH1F("tpcWXInSurveyXZ","tpcWXInSurveyXZ",201,-1.005,1.005);
121  tpcWXInSurveyXZ->SetXTitle("Angle (mrad)");
122 
123  TH1F *tpcWYInSurveyXY = new TH1F("tpcWYInSurveyXY","tpcWYInSurveyXY",201,-1.005,1.005);
124  tpcWYInSurveyXY->SetXTitle("Angle (mrad)");
125  TH1F *tpcWYInSurveyYZ = new TH1F("tpcWYInSurveyYZ","tpcWYInSurveyYZ",201,-1.005,1.005);
126  tpcWYInSurveyYZ->SetXTitle("Angle (mrad)");
127 
128  TH1F *tpcWZInSurveyXZ = new TH1F("tpcWZInSurveyXZ","tpcWZInSurveyXZ",201,-1.005,1.005);
129  TH1F *tpcWZInSurveyYZ = new TH1F("tpcWZInSurveyYZ","tpcWZInSurveyYZ",201,-1.005,1.005);
130 
131  TH1F *tpcEXInSurveyXY = new TH1F("tpcEXInSurveyXY","tpcEXInSurveyXY",201,-1.005,1.005);
132  TH1F *tpcEXInSurveyXZ = new TH1F("tpcEXInSurveyXZ","tpcEXInSurveyXZ",201,-1.005,1.005);
133 
134  TH1F *tpcEYInSurveyXY = new TH1F("tpcEYInSurveyXY","tpcEYInSurveyXY",201,-1.005,1.005);
135  tpcEYInSurveyXY->SetXTitle("Angle (mrad)");
136  TH1F *tpcEYInSurveyYZ = new TH1F("tpcEYInSurveyYZ","tpcEYInSurveyYZ",201,-1.005,1.005);
137  tpcEYInSurveyYZ->SetXTitle("Angle (mrad)");
138 
139  TH1F *tpcEZInSurveyXZ = new TH1F("tpcEZInSurveyXZ","tpcEZInSurveyXZ",201,-1.005,1.005);
140  TH1F *tpcEZInSurveyYZ = new TH1F("tpcEZInSurveyYZ","tpcEZInSurveyYZ",201,-1.005,1.005);
141 
142 
143 
144  //
145  // To calculate the center of the TPC (perfect cylinder, 12 perfectly symetric sectors at each end)
146  // Calculate the midle of a line going from one end to the other through the center of the TPC
147  // W1-E17, W2-E16, W3-E15, W4-E14, W5-E13
148  // W6-E24, W7-E23, W8-E22, W9-E21, W10-E20, W11-E19, W12-E18
149  for(int iSec=2; iSec<=5; iSec++){ // No WAO[1] data
150  tempVect.SetXYZ(0.,0.,0.);
151  tempVect += WAO[iSec]; tempVect += EAO[18-iSec]; tempVect *= 0.5; tempVect *= 1000;
152  tpcCenterInSurvey->Fill(tempVect.X(),tempVect.Y(),tempVect.Z());
153  tpcCenterInSurveyX->Fill(tempVect.X());
154  tpcCenterInSurveyY->Fill(tempVect.Y());
155  tpcCenterInSurveyZ->Fill(tempVect.Z());
156  printf("%f\n",tempVect.X());
157  }
158  for(int iSec=6; iSec<=12; iSec++){
159  tempVect.SetXYZ(0.,0.,0.);
160  tempVect += WAO[iSec]; tempVect += EAO[30-iSec]; tempVect *= 0.5; tempVect *= 1000;
161  tpcCenterInSurvey->Fill(tempVect.X(),tempVect.Y(),tempVect.Z());
162  tpcCenterInSurveyX->Fill(tempVect.X());
163  tpcCenterInSurveyY->Fill(tempVect.Y());
164  tpcCenterInSurveyZ->Fill(tempVect.Z());
165  printf("%f\n",tempVect.X());
166  }
167 
168  for(int iSec=1; iSec<=5; iSec++){
169  tempVect.SetXYZ(0.,0.,0.);
170  tempVect += WBO[iSec]; tempVect += EBO[18-iSec]; tempVect *= 0.5; tempVect *= 1000;
171  tpcCenterInSurvey->Fill(tempVect.X(),tempVect.Y(),tempVect.Z());
172  tpcCenterInSurveyX->Fill(tempVect.X());
173  tpcCenterInSurveyY->Fill(tempVect.Y());
174  tpcCenterInSurveyZ->Fill(tempVect.Z());
175  printf("%f\n",tempVect.X());
176  }
177  for(int iSec=6; iSec<=12; iSec++){
178  if(iSec!=11){ // No EBO[19] data
179  tempVect.SetXYZ(0.,0.,0.);
180  tempVect += WBO[iSec]; tempVect += EBO[30-iSec]; tempVect *= 0.5; tempVect *= 1000;
181  tpcCenterInSurvey->Fill(tempVect.X(),tempVect.Y(),tempVect.Z());
182  tpcCenterInSurveyX->Fill(tempVect.X());
183  tpcCenterInSurveyY->Fill(tempVect.Y());
184  tpcCenterInSurveyZ->Fill(tempVect.Z());
185  printf("%f\n",tempVect.X());
186  }
187  }
188 
189  for(int iSec=1; iSec<=5; iSec++){
190  if(iSec!=3&&iSec!=4){ // No EDO[19] data
191  tempVect.SetXYZ(0.,0.,0.);
192  tempVect += WDO[iSec]; tempVect += EDO[18-iSec]; tempVect *= 0.5; tempVect *= 1000;
193  tpcCenterInSurvey->Fill(tempVect.X(),tempVect.Y(),tempVect.Z());
194  tpcCenterInSurveyX->Fill(tempVect.X());
195  tpcCenterInSurveyY->Fill(tempVect.Y());
196  tpcCenterInSurveyZ->Fill(tempVect.Z());
197  printf("%f\n",tempVect.X());
198  }
199  }
200  for(int iSec=6; iSec<=12; iSec++){
201  if(iSec!=8&&iSec!=9&&iSec!=11&&iSec!=12){ // No EDO[19] data
202  tempVect.SetXYZ(0.,0.,0.);
203  tempVect += WDO[iSec]; tempVect += EDO[30-iSec]; tempVect *= 0.5; tempVect *= 1000;
204  tpcCenterInSurvey->Fill(tempVect.X(),tempVect.Y(),tempVect.Z());
205  tpcCenterInSurveyX->Fill(tempVect.X());
206  tpcCenterInSurveyY->Fill(tempVect.Y());
207  tpcCenterInSurveyZ->Fill(tempVect.Z());
208  printf("%f\n",tempVect.X());
209  }
210  }
211 
212  //
213  // To calculate the angle TPC with respect to the Surveynet,
214  // Calculate the normal vector of a line in the endcap plane
215  //
216  // A targets
217  for(int iSec=1; iSec<=10; iSec++){
218  if(WAO[iSec]!=nullVect){
219  for(int jSec=iSec+1; jSec<=11; jSec++){
220  if(WAO[jSec]!=nullVect){
221  tempVect.SetXYZ(0.,0.,0.);
222  tempVect += WAO[iSec]; tempVect -= WAO[jSec];
223  for(int kSec=jSec+1; kSec<=12; kSec++){
224  if(WAO[kSec]!=nullVect){
225  crossVect.SetXYZ(0.,0.,0.);
226  crossVect += WAO[jSec]; crossVect -= WAO[kSec];
227  crossVect = crossVect.Cross(tempVect);
228  uniVect = crossVect.Unit();
229  //printf("(%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,uniVect.X(),uniVect.Y(),uniVect.Z());
230 
231  tAngle = TMath::ATan2(uniVect.X(),uniVect.Y());
232  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
233  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
234  tAngle *= 1e3;
235  if (tAngle>10000) printf("WAO (%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,uniVect.X(),uniVect.Y(),uniVect.Z());
236  //printf("In XY plane, Angle to Y %f\n",tAngle);
237  tpcWYInSurveyXY->Fill(tAngle);
238  tAngle = TMath::ATan2(uniVect.Z(),uniVect.Y());
239  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
240  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
241  tAngle *= 1e3;
242  if (tAngle>10000) printf("WAO (%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,uniVect.X(),uniVect.Y(),uniVect.Z());
243  //printf("In YZ plane, Angle to Y %f\n",tAngle);
244  tpcWYInSurveyYZ->Fill(tAngle);
245  }
246  }
247  }
248  }
249  }
250  }
251 
252  for(int iSec=1; iSec<=9; iSec++){
253  if(WAO[iSec]!=nullVect){
254  for(int jSec=iSec+1; jSec<=10; jSec++){
255  if(WAO[jSec]!=nullVect){
256  tempVect.SetXYZ(0.,0.,0.);
257  tempVect += WAO[iSec]; tempVect -= WAO[jSec];
258  for(int kSec=jSec+1; kSec<=11; kSec++){
259  if(WAO[kSec]!=nullVect){
260  for(int lSec=kSec+1; lSec<=12; lSec++){
261  if(WAO[lSec]!=nullVect){
262  crossVect.SetXYZ(0.,0.,0.);
263  crossVect += WAO[kSec]; crossVect -= WAO[lSec];
264  crossVect = crossVect.Cross(tempVect);
265  uniVect = crossVect.Unit();
266 
267  tAngle = TMath::ATan2(uniVect.X(),uniVect.Y());
268  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
269  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
270  tAngle *= 1e3;
271  if (tAngle>10000) printf("WAO (%d,%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,lSec,uniVect.X(),uniVect.Y(),uniVect.Z());
272  //printf("In XY plane, Angle to Y %f\n",tAngle);
273  tpcWYInSurveyXY->Fill(tAngle);
274  tAngle = TMath::ATan2(uniVect.Z(),uniVect.Y());
275  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
276  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
277  tAngle *= 1e3;
278  if (tAngle>10000) printf("WAO (%d,%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,lSec,uniVect.X(),uniVect.Y(),uniVect.Z());
279  //printf("In YZ plane, Angle to Y %f\n",tAngle);
280  tpcWYInSurveyYZ->Fill(tAngle);
281  }
282  }
283  }
284  }
285  }
286  }
287  }
288  }
289 
290  // B targets
291  for(int iSec=1; iSec<=10; iSec++){
292  if(WBO[iSec]!=nullVect){
293  for(int jSec=iSec+1; jSec<=11; jSec++){
294  if(WBO[jSec]!=nullVect){
295  tempVect.SetXYZ(0.,0.,0.);
296  tempVect += WBO[iSec]; tempVect -= WBO[jSec];
297  for(int kSec=jSec+1; kSec<=12; kSec++){
298  if(WBO[kSec]!=nullVect){
299  crossVect.SetXYZ(0.,0.,0.);
300  crossVect += WBO[jSec]; crossVect -= WBO[kSec];
301  crossVect = crossVect.Cross(tempVect);
302  uniVect = crossVect.Unit();
303 
304  tAngle = TMath::ATan2(uniVect.X(),uniVect.Y());
305  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
306  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
307  tAngle *= 1e3;
308  if (tAngle>10000) printf("WBO (%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,uniVect.X(),uniVect.Y(),uniVect.Z());
309  //printf("In XY plane, Angle to Y %f\n",tAngle);
310  tpcWYInSurveyXY->Fill(tAngle);
311  tAngle = TMath::ATan2(uniVect.Z(),uniVect.Y());
312  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
313  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
314  tAngle *= 1e3;
315  if (tAngle>10000) printf("WBO (%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,uniVect.X(),uniVect.Y(),uniVect.Z());
316  //printf("In YZ plane, Angle to Y %f\n",tAngle);
317  tpcWYInSurveyYZ->Fill(tAngle);
318  }
319  }
320  }
321  }
322  }
323  }
324 
325  for(int iSec=1; iSec<=9; iSec++){
326  if(WBO[iSec]!=nullVect){
327  for(int jSec=iSec+1; jSec<=10; jSec++){
328  if(WBO[jSec]!=nullVect){
329  tempVect.SetXYZ(0.,0.,0.);
330  tempVect += WBO[iSec]; tempVect -= WBO[jSec];
331  for(int kSec=jSec+1; kSec<=11; kSec++){
332  if(WBO[kSec]!=nullVect){
333  for(int lSec=kSec+1; lSec<=12; lSec++){
334  if(WBO[lSec]!=nullVect){
335  crossVect.SetXYZ(0.,0.,0.);
336  crossVect += WBO[kSec]; crossVect -= WBO[lSec];
337  crossVect = crossVect.Cross(tempVect);
338  uniVect = crossVect.Unit();
339 
340  tAngle = TMath::ATan2(uniVect.X(),uniVect.Y());
341  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
342  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
343  tAngle *= 1e3;
344  if (tAngle>10000) printf("WBO (%d,%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,lSec,uniVect.X(),uniVect.Y(),uniVect.Z());
345  //printf("In XY plane, Angle to Y %f\n",tAngle);
346  tpcWYInSurveyXY->Fill(tAngle);
347  tAngle = TMath::ATan2(uniVect.Z(),uniVect.Y());
348  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
349  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
350  tAngle *= 1e3;
351  if (tAngle>10000) printf("WBO (%d,%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,lSec,uniVect.X(),uniVect.Y(),uniVect.Z());
352  //printf("In YZ plane, Angle to Y %f\n",tAngle);
353  tpcWYInSurveyYZ->Fill(tAngle);
354  }
355  }
356  }
357  }
358  }
359  }
360  }
361  }
362 
363  // D targets
364  for(int iSec=1; iSec<=10; iSec++){
365  if(WDO[iSec]!=nullVect){
366  for(int jSec=iSec+1; jSec<=11; jSec++){
367  if(WDO[jSec]!=nullVect){
368  tempVect.SetXYZ(0.,0.,0.);
369  tempVect += WDO[iSec]; tempVect -= WDO[jSec];
370  for(int kSec=jSec+1; kSec<=12; kSec++){
371  if(WDO[kSec]!=nullVect){
372  crossVect.SetXYZ(0.,0.,0.);
373  crossVect += WDO[jSec]; crossVect -= WDO[kSec];
374  crossVect = crossVect.Cross(tempVect);
375  uniVect = crossVect.Unit();
376 
377  tAngle = TMath::ATan2(uniVect.X(),uniVect.Y());
378  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
379  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
380  tAngle *= 1e3;
381  if (tAngle>10000) printf("WDO (%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,uniVect.X(),uniVect.Y(),uniVect.Z());
382  //printf("In XY plane, Angle to Y %f\n",tAngle);
383  tpcWYInSurveyXY->Fill(tAngle);
384  tAngle = TMath::ATan2(uniVect.Z(),uniVect.Y());
385  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
386  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
387  tAngle *= 1e3;
388  if (tAngle>10000) printf("WDO (%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,uniVect.X(),uniVect.Y(),uniVect.Z());
389  //printf("In YZ plane, Angle to Y %f\n",tAngle);
390  tpcWYInSurveyYZ->Fill(tAngle);
391  }
392  }
393  }
394  }
395  }
396  }
397 
398  for(int iSec=1; iSec<=9; iSec++){
399  if(WDO[iSec]!=nullVect){
400  for(int jSec=iSec+1; jSec<=10; jSec++){
401  if(WDO[jSec]!=nullVect){
402  tempVect.SetXYZ(0.,0.,0.);
403  tempVect += WDO[iSec]; tempVect -= WDO[jSec];
404  for(int kSec=jSec+1; kSec<=11; kSec++){
405  if(WDO[kSec]!=nullVect){
406  for(int lSec=kSec+1; lSec<=12; lSec++){
407  if(WDO[lSec]!=nullVect){
408  crossVect.SetXYZ(0.,0.,0.);
409  crossVect += WDO[kSec]; crossVect -= WDO[lSec];
410  crossVect = crossVect.Cross(tempVect);
411  uniVect = crossVect.Unit();
412 
413  tAngle = TMath::ATan2(uniVect.X(),uniVect.Y());
414  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
415  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
416  tAngle *= 1e3;
417  if (tAngle>10000) printf("WDO (%d,%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,lSec,uniVect.X(),uniVect.Y(),uniVect.Z());
418  //printf("In XY plane, Angle to Y %f\n",tAngle);
419  tpcWYInSurveyXY->Fill(tAngle);
420  tAngle = TMath::ATan2(uniVect.Z(),uniVect.Y());
421  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
422  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
423  tAngle *= 1e3;
424  if (tAngle>10000) printf("WDO (%d,%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,lSec,uniVect.X(),uniVect.Y(),uniVect.Z());
425  //printf("In YZ plane, Angle to Y %f\n",tAngle);
426  tpcWYInSurveyYZ->Fill(tAngle);
427  }
428  }
429  }
430  }
431  }
432  }
433  }
434  }
435 
436 
437  // East
438  // A targets
439  for(int iSec=13; iSec<=22; iSec++){
440  if(EAO[iSec]!=nullVect){
441  for(int jSec=iSec+1; jSec<=23; jSec++){
442  if(EAO[jSec]!=nullVect){
443  tempVect.SetXYZ(0.,0.,0.);
444  tempVect += EAO[iSec]; tempVect -= EAO[jSec];
445  for(int kSec=jSec+1; kSec<=24; kSec++){
446  if(EAO[kSec]!=nullVect){
447  crossVect.SetXYZ(0.,0.,0.);
448  crossVect += EAO[jSec]; crossVect -= EAO[kSec];
449  crossVect = crossVect.Cross(tempVect);
450  uniVect = crossVect.Unit();
451 
452  tAngle = TMath::ATan2(uniVect.X(),uniVect.Y());
453  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
454  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
455  tAngle *= 1e3;
456  if (tAngle>10000) printf("EAO (%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,uniVect.X(),uniVect.Y(),uniVect.Z());
457  //printf("In XY plane, Angle to Y %f\n",tAngle);
458  tpcEYInSurveyXY->Fill(tAngle);
459  tAngle = TMath::ATan2(uniVect.Z(),uniVect.Y());
460  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
461  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
462  tAngle *= 1e3;
463  if (tAngle>10000) printf("EAO (%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,uniVect.X(),uniVect.Y(),uniVect.Z());
464  //printf("In YZ plane, Angle to Y %f\n",tAngle);
465  tpcEYInSurveyYZ->Fill(tAngle);
466  }
467  }
468  }
469  }
470  }
471  }
472 
473  for(int iSec=13; iSec<=21; iSec++){
474  if(EAO[iSec]!=nullVect){
475  for(int jSec=iSec+1; jSec<=22; jSec++){
476  if(EAO[jSec]!=nullVect){
477  tempVect.SetXYZ(0.,0.,0.);
478  tempVect += EAO[iSec]; tempVect -= EAO[jSec];
479  for(int kSec=jSec+1; kSec<=23; kSec++){
480  if(EAO[kSec]!=nullVect){
481  for(int lSec=kSec+1; lSec<=24; lSec++){
482  if(EAO[lSec]!=nullVect){
483  crossVect.SetXYZ(0.,0.,0.);
484  crossVect += EAO[kSec]; crossVect -= EAO[lSec];
485  crossVect = crossVect.Cross(tempVect);
486  uniVect = crossVect.Unit();
487 
488  tAngle = TMath::ATan2(uniVect.X(),uniVect.Y());
489  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
490  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
491  tAngle *= 1e3;
492  if (tAngle>10000) printf("EAO (%d,%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,lSec,uniVect.X(),uniVect.Y(),uniVect.Z());
493  //printf("In XY plane, Angle to Y %f\n",tAngle);
494  tpcEYInSurveyXY->Fill(tAngle);
495  tAngle = TMath::ATan2(uniVect.Z(),uniVect.Y());
496  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
497  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
498  tAngle *= 1e3;
499  if (tAngle>10000) printf("EAO (%d,%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,lSec,uniVect.X(),uniVect.Y(),uniVect.Z());
500  //printf("In YZ plane, Angle to Y %f\n",tAngle);
501  tpcEYInSurveyYZ->Fill(tAngle);
502  }
503  }
504  }
505  }
506  }
507  }
508  }
509  }
510 
511  // B targets
512  for(int iSec=13; iSec<=22; iSec++){
513  if(EBO[iSec]!=nullVect){
514  for(int jSec=iSec+1; jSec<=23; jSec++){
515  if(EBO[jSec]!=nullVect){
516  tempVect.SetXYZ(0.,0.,0.);
517  tempVect += EBO[iSec]; tempVect -= EBO[jSec];
518  for(int kSec=jSec+1; kSec<=24; kSec++){
519  if(EBO[kSec]!=nullVect){
520  crossVect.SetXYZ(0.,0.,0.);
521  crossVect += EBO[jSec]; crossVect -= EBO[kSec];
522  crossVect = crossVect.Cross(tempVect);
523  uniVect = crossVect.Unit();
524 
525  tAngle = TMath::ATan2(uniVect.X(),uniVect.Y());
526  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
527  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
528  tAngle *= 1e3;
529  if (tAngle>10000) printf("EBO (%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,uniVect.X(),uniVect.Y(),uniVect.Z());
530  //printf("In XY plane, Angle to Y %f\n",tAngle);
531  tpcEYInSurveyXY->Fill(tAngle);
532  tAngle = TMath::ATan2(uniVect.Z(),uniVect.Y());
533  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
534  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
535  tAngle *= 1e3;
536  if (tAngle>10000) printf("EBO (%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,uniVect.X(),uniVect.Y(),uniVect.Z());
537  //printf("In YZ plane, Angle to Y %f\n",tAngle);
538  tpcEYInSurveyYZ->Fill(tAngle);
539  }
540  }
541  }
542  }
543  }
544  }
545 
546  for(int iSec=13; iSec<=21; iSec++){
547  if(EBO[iSec]!=nullVect){
548  for(int jSec=iSec+1; jSec<=22; jSec++){
549  if(EBO[jSec]!=nullVect){
550  tempVect.SetXYZ(0.,0.,0.);
551  tempVect += EBO[iSec]; tempVect -= EBO[jSec];
552  for(int kSec=jSec+1; kSec<=23; kSec++){
553  if(EBO[kSec]!=nullVect){
554  for(int lSec=kSec+1; lSec<=24; lSec++){
555  if(EBO[lSec]!=nullVect){
556  crossVect.SetXYZ(0.,0.,0.);
557  crossVect += EBO[kSec]; crossVect -= EBO[lSec];
558  crossVect = crossVect.Cross(tempVect);
559  uniVect = crossVect.Unit();
560 
561  tAngle = TMath::ATan2(uniVect.X(),uniVect.Y());
562  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
563  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
564  tAngle *= 1e3;
565  if (tAngle>10000) printf("EBO (%d,%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,lSec,uniVect.X(),uniVect.Y(),uniVect.Z());
566  //printf("In XY plane, Angle to Y %f\n",tAngle);
567  tpcEYInSurveyXY->Fill(tAngle);
568  tAngle = TMath::ATan2(uniVect.Z(),uniVect.Y());
569  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
570  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
571  tAngle *= 1e3;
572  if (tAngle>10000) printf("EBO (%d,%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,lSec,uniVect.X(),uniVect.Y(),uniVect.Z());
573  //printf("In YZ plane, Angle to Y %f\n",tAngle);
574  tpcEYInSurveyYZ->Fill(tAngle);
575  }
576  }
577  }
578  }
579  }
580  }
581  }
582  }
583 
584  // D targets
585  for(int iSec=13; iSec<=22; iSec++){
586  if(EDO[iSec]!=nullVect){
587  for(int jSec=iSec+1; jSec<=23; jSec++){
588  if(EDO[jSec]!=nullVect){
589  tempVect.SetXYZ(0.,0.,0.);
590  tempVect += EDO[iSec]; tempVect -= EDO[jSec];
591  for(int kSec=jSec+1; kSec<=24; kSec++){
592  if(EDO[kSec]!=nullVect){
593  crossVect.SetXYZ(0.,0.,0.);
594  crossVect += EDO[jSec]; crossVect -= EDO[kSec];
595  crossVect = crossVect.Cross(tempVect);
596  uniVect = crossVect.Unit();
597 
598  tAngle = TMath::ATan2(uniVect.X(),uniVect.Y());
599  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
600  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
601  tAngle *= 1e3;
602  if (tAngle>10000) printf("EDO (%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,uniVect.X(),uniVect.Y(),uniVect.Z());
603  //printf("In XY plane, Angle to Y %f\n",tAngle);
604  tpcEYInSurveyXY->Fill(tAngle);
605  tAngle = TMath::ATan2(uniVect.Z(),uniVect.Y());
606  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
607  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
608  tAngle *= 1e3;
609  if (tAngle>10000) printf("EDO (%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,uniVect.X(),uniVect.Y(),uniVect.Z());
610  //printf("In YZ plane, Angle to Y %f\n",tAngle);
611  tpcEYInSurveyYZ->Fill(tAngle);
612  }
613  }
614  }
615  }
616  }
617  }
618 
619  for(int iSec=13; iSec<=21; iSec++){
620  if(EDO[iSec]!=nullVect){
621  for(int jSec=iSec+1; jSec<=22; jSec++){
622  if(EDO[jSec]!=nullVect){
623  tempVect.SetXYZ(0.,0.,0.);
624  tempVect += EDO[iSec]; tempVect -= EDO[jSec];
625  for(int kSec=jSec+1; kSec<=23; kSec++){
626  if(EDO[kSec]!=nullVect){
627  for(int lSec=kSec+1; lSec<=24; lSec++){
628  if(EDO[lSec]!=nullVect){
629  crossVect.SetXYZ(0.,0.,0.);
630  crossVect += EDO[kSec]; crossVect -= EDO[lSec];
631  crossVect = crossVect.Cross(tempVect);
632  uniVect = crossVect.Unit();
633 
634  tAngle = TMath::ATan2(uniVect.X(),uniVect.Y());
635  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
636  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
637  tAngle *= 1e3;
638  if (tAngle>10000) printf("EDO (%d,%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,lSec,uniVect.X(),uniVect.Y(),uniVect.Z());
639  //printf("In XY plane, Angle to Y %f\n",tAngle);
640  tpcEYInSurveyXY->Fill(tAngle);
641  tAngle = TMath::ATan2(uniVect.Z(),uniVect.Y());
642  if (tAngle>TMath::Pi()/2.) tAngle -= TMath::Pi();
643  if (tAngle<-TMath::Pi()/2.) tAngle += TMath::Pi();
644  tAngle *= 1e3;
645  if (tAngle>10000) printf("EDO (%d,%d,%d,%d) (%f,%f,%f)\n",iSec,jSec,kSec,lSec,uniVect.X(),uniVect.Y(),uniVect.Z());
646  //printf("In YZ plane, Angle to Y %f\n",tAngle);
647  tpcEYInSurveyYZ->Fill(tAngle);
648  }
649  }
650  }
651  }
652  }
653  }
654  }
655  }
656 
657 
658 
659 /* for(int iSec=13; iSec<=15; iSec++){ */
660 /* tempVect.SetXYZ(0.,0.,0.); */
661 /* tempVect += EAO[iSec]; tempVect -= EAO[iSec+6]; */
662 /* crossVect.SetXYZ(0.,0.,0.); */
663 /* crossVect += EAO[iSec+3]; crossVect -= EAO[iSec+3+6]; */
664 /* crossVect = crossVect.Cross(tempVect); */
665 
666 /* uniVect = crossVect.Unit(); */
667 /* uniVect.SetX(0.); */
668 /* tAngle = uniVect.Angle(axisSurveyY); */
669 /* tAngle *= 1e3; */
670 /* printf("In YZ plane, Angle to Y %f\n",tAngle); */
671 /* uniVect = crossVect.Unit(); */
672 /* uniVect.SetZ(0.); */
673 /* tAngle = uniVect.Angle(axisSurveyY); */
674 /* tAngle *= 1e3; */
675 /* printf("In XY plane, Angle to Y %f\n",tAngle); */
676 
677 /* } */
678 
679 }
Generated by   doxygen 1.8.5