StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
CheckLasers.C
1 #include "LaserBeams.h"
2 #include "LASERINO.h"
3 void CheckLasers() {
4  cout << "NoLocals " << NoLocals << "\tNoBeams " << NoBeams << endl;
5  const Int_t numSectors = 24;
6  Double_t beta = 0;
7  Double_t dBeta = 720./numSectors;
8  Int_t sector = 0;
9  Int_t half = 0;
10  TGeoHMatrix TpcHalf[2];
11  Double_t rotHalfs[18] = {
12  0, -1, 0, -1, 0, 0, 0, 0, -1,// sector 1-12
13  0, 1, 0, -1, 0, 0, 0, 0, 1 // sector 13-24
14  };
15  for (half = 0; half <2; half++) TpcHalf[half].SetRotation(&rotHalfs[9*half]);
16  TGeoHMatrix RotSec[24];
17  for (sector = 1; sector <= numSectors; sector++) {
18  if (sector > 12) beta = (numSectors-sector)*dBeta;
19  else beta = sector*dBeta;
20  RotSec[sector-1].RotateZ(-beta);
21  }
22 
23  for (Int_t raft = 1; raft <= 14; raft++) {
24  // for (Int_t raft = 9; raft <= 9; raft++) {
25  TRVector mX;
26  TRMatrix A(0,6);
27  for (Int_t bundle = 1; bundle <= 6; bundle++) {
28  for (Int_t mirror = 1; mirror <= 7; mirror++) {
29  Int_t RBM = 100*raft + 10*bundle + mirror;
30  if (RBM == 121 || // +1
31  RBM == 122 || // +1
32  RBM == 124 || // +1
33  RBM == 126 || // +1
34  RBM == 131 || // +1
35  RBM == 132 ) // +1
36  continue;
37  if (RBM == 211 || // +1
38  RBM == 221 || // +1
39  RBM == 225 || // +1
40  RBM == 231 || // +1
41  RBM == 244 || // +1
42  RBM == 246 || // +1
43  RBM == 265) // +1
44  continue;
45  if (RBM == 311 || // +1
46  RBM == 331 || // +1
47 // RBM == 335 || // -1
48  RBM == 337 ) // +1
49  continue;
50  if (raft == 5) continue;
51  Int_t l = -1;
52  for (Int_t local = 0; local < NoLocals; local++) {
53  if (Locals[local].Raft == raft && Locals[local].Bundle == bundle && Locals[local].Mirror == mirror) {l = local; break;}
54  }
55  Int_t k = -1;
56  for (Int_t global = 0; global < NoBeams; global++) {
57  if (LaserBeams[global].Raft == raft && LaserBeams[global].Bundle == bundle && LaserBeams[global].Mirror == mirror) {k = global; break;}
58  }
59  if (l < 0 || k < 0) continue;
60  // cout << "Found match for sector " << LaserBeams[k].Sector << " raft " << raft << " bundle " << bundle << " mirror " << mirror << endl;
61  Double_t *xg = &LaserBeams[k].X;// TRVector XG(3,xg); cout << " XG " << XG << endl;
62  Double_t xG[3];
63  sector = LaserBeams[k].Sector;
64  half = 0;
65  if (sector > 12) half = 1;
66  TGeoHMatrix RotM = RotSec[sector-1] * TpcHalf[half];
67  // TGeoHMatrix RotM = TpcHalf[half] * RotSec[sector-1];
68  // cout << " Raft " << raft << " Sector " << sector << " Bundle " << bundle << " Mirror " << mirror << endl;
69  // RotM.Print(); //TRVector xgT(3,xg); cout << "xg " << xgT;
70  RotM.MasterToLocal(xg,xG); //TRVector xGT(3,xG); cout << "xG " << xGT;
71  Double_t *xl = &Locals[l].X;
72  Double_t xL[3] = {0.1*xl[0], 0.1*xl[1], 0.1*xl[2]};// TRVector XL(3,xL); cout << "XL " << XL;
73  /*
74  ( 1 -gamma beta )
75  DRT = ( gamma 1 -alpha);
76  ( -beta alpha 1 )
77  */
78  Double_t a[18] = {
79  0., xL[2], -xL[1], 1., 0., 0.,
80  -xL[2], 0., xL[0], 0., 1., 0.,
81  xL[1], -xL[0], 0., 0., 0., 1.
82  };
83  for (Int_t i = 0; i < 3; i++) {
84  Double_t x = xG[i] - xL[i];
85  mX.AddRow(&x);
86  A.AddRow(&a[6*i]);
87  }
88  }
89  }
90  // cout << "X "; mX.Print();
91  // cout << "A "; A.Print();
92  if (A.GetNrows() < 1) {
93  cout << "Raft " << raft << " is empty" << endl;
94  continue;
95  }
96  cout << "Raft " << raft << " Sector = " << sector << " =========================" << endl;
97  TRVector AmX(A,TRArray::kATxB,mX); // cout << "AmX\t" << AmX << endl;
98  TRSymMatrix S(A,TRArray::kATxA); // cout << "S\t" << S << endl;
99  TRSymMatrix SInv(S,TRArray::kInverted); // cout << "SInv " << SInv << endl;
100  TRVector X(SInv,TRArray::kSxA,AmX); // cout << "X " << X << endl;
101  TGeoHMatrix T(Form("Raft%0i",raft));
102  // normolize
103  Double_t *xx = X.GetArray();
104  T.RotateX(xx[0]*180/TMath::Pi());
105  T.RotateY(xx[1]*180/TMath::Pi());
106  T.RotateZ(xx[2]*180/TMath::Pi());
107  T.SetTranslation(&xx[3]);
108  cout << Form("{%i,%i",raft,sector);
109  for (Int_t i = 0; i < 6; i++) cout << Form(",%15.5f",xx[i]);
110  cout << "}," << endl;
111  // cout << "Determinant " << T.Determinant() << endl;
112  T.Print();
113  TGeoHMatrix RotM = RotSec[sector-1] * TpcHalf[half] * T;
114  // TGeoHMatrix RotM = TpcHalf[half] * RotSec[sector-1] * T;
115  RotM.Print();
116  // Check table
117  for (Int_t bundle = 1; bundle <= 6; bundle++) {
118  for (Int_t mirror = 1; mirror <= 7; mirror++) {
119  Int_t l = -1;
120  for (Int_t local = 0; local < NoLocals; local++) {
121  if (Locals[local].Raft == raft && Locals[local].Bundle == bundle && Locals[local].Mirror == mirror) {l = local; break;}
122  }
123  Int_t k = -1;
124  for (Int_t global = 0; global < NoBeams; global++) {
125  if (LaserBeams[global].Raft == raft && LaserBeams[global].Bundle == bundle && LaserBeams[global].Mirror == mirror) {k = global; break;}
126  }
127  if (l < 0 || k < 0) continue;
128  sector = LaserBeams[k].Sector;
129  // cout << "Found match for raft " << raft << " bundle " << bundle << " mirror " << mirror << endl;
130  Double_t *xg = &LaserBeams[k].X; // TRVector XG(3,xg); cout << " XG " << XG;
131  Double_t *xl = &Locals[l].X;
132  Double_t xL[3] = {0.1*xl[0],0.1*xl[1],0.1*xl[2]}; // TRVector XL(3,xL); cout << " xL " << XL;
133  Double_t xG[3];
134 
135  RotM.LocalToMaster(xL,xG); // TRVector XGG(3,xG); cout << " xG " << XGG;
136  Double_t dev = 0;
137  for (Int_t i = 0; i < 3; i++) {
138  Double_t dif = (xG[i] - xg[i]);
139  dev += dif*dif;
140  }
141  dev = TMath::Sqrt(dev);
142  if (dev > 0.0000) {
143  cout << Form("S%02iR%02iB%iM%i dev = %f",sector,raft,bundle,mirror,dev);
144  cout << Form(" dX/dY/dZ = %f/%f/%f",xG[0] - xg[0],xG[1] - xg[1],xG[2] - xg[2]);
145  if (dev > 0.005) cout << " ===========";
146  cout << endl;
147  }
148  }
149  }
150  }
151 }
Generated by   doxygen 1.8.5