FairRoot/PandaRoot
PndRichCalDb.cxx
Go to the documentation of this file.
1 // -------------------------------------------------------------------------
2 // ----- PndRichCalDb source file -----
3 // ----- Created 12/01/04 by M. Al/Turany (FairField.cxx) -----
4 // -------------------------------------------------------------------------
5 
6 
7 #include <iomanip>
8 #include <iostream>
9 #include <fstream>
10 #include "stdlib.h"
11 
12 #include "TArrayF.h"
13 #include "TFile.h"
14 #include "TMath.h"
15 #include "TObjArray.h"
16 
17 #include "PndRichCalDb.h"
18 #include "PndRichCalDbData.h"
19 #include "PndRichCalDbPar.h"
20 
21 using namespace std;
22 
23 // ------------- Default constructor ----------------------------------
25  : FairField(),
26  fFileName(""),
27  fPmin(0), fPmax(0), fPstep(0),
28  fXmin(0), fXmax(0), fXstep(0),
29  fYmin(0), fYmax(0), fYstep(0),
30  fTmin(0), fTmax(0), fTstep(0),
31  fFmin(0), fFmax(0), fFstep(0),
32  fNp(0), fNx(0), fNy(0), fNt(0), fNf(0),
33  fBetaMean(NULL), fBetaSig(NULL), fBetaEff(NULL)
34 {
35  SetName("");
36  fType = 0;
37  fPmin=1; fPmax=11; fPstep=10;
38  fXmin=2; fXmax=12; fXstep=10;
39  fYmin=3; fYmax=13; fYstep=10;
40  fTmin=4; fTmax=14; fTstep=10;
41  fFmin=5; fFmax=15; fFstep=10;
42  fNp=2; fNx=2; fNy=2; fNt=2; fNf=2;
43  Int_t size = fNp*fNx*fNy*fNt*fNf;
44  fBetaSig = new TArrayF(size);
45  fBetaMean = new TArrayF(size);
46  fBetaEff = new TArrayF(size);
47 
48 }
49 // ------------------------------------------------------------------------
50 
51 
52 
53 // ------------- Standard constructor ---------------------------------
54 PndRichCalDb::PndRichCalDb(const char* mapName, const char* fileType)
55  : FairField(mapName),
56  fFileName(""),
57  fPmin(0), fPmax(0), fPstep(0),
58  fXmin(0), fXmax(0), fXstep(0),
59  fYmin(0), fYmax(0), fYstep(0),
60  fTmin(0), fTmax(0), fTstep(0),
61  fFmin(0), fFmax(0), fFstep(0),
62  fNp(0), fNx(0), fNy(0), fNt(0), fNf(0),
63  fBetaMean(NULL), fBetaSig(NULL), fBetaEff(NULL)
64 {
65  SetName(mapName);
66  TString dir = getenv("VMCWORKDIR");
67  fFileName = dir + "/detectors/rich/" + mapName;
68  if ( fileType[0] == 'R' ) fFileName += ".root";
69  else fFileName += ".dat";
70  fType = 1;
71 }
72 // ------------------------------------------------------------------------
73 
74 
75 
76 // ------------ Constructor from PndRichCalDbPar --------------------------
78  : FairField(),
79  fFileName(""),
80  fPmin(0), fPmax(0), fPstep(0),
81  fXmin(0), fXmax(0), fXstep(0),
82  fYmin(0), fYmax(0), fYstep(0),
83  fTmin(0), fTmax(0), fTstep(0),
84  fFmin(0), fFmax(0), fFstep(0),
85  fNp(0), fNx(0), fNy(0), fNt(0), fNf(0),
86  fBetaMean(NULL), fBetaSig(NULL), fBetaEff(NULL)
87 {
88  fType = 1;
89  if ( ! calDbPar ) {
90  LOG(ERROR) << "PndRichCalDb::PndRichCalDb: empty parameter container!";
91  SetName("");
92  fType = -1;
93  }
94  else {
95  TString Name=GetName();
96  calDbPar->MapName(Name);
97  TString dir = getenv("VMCWORKDIR");
98  fFileName = dir + "/detectors/rich/" + Name + ".root";
99  fFileName = Name + ".root";
100  fType = calDbPar->GetType();
101  }
102 }
103 // ------------------------------------------------------------------------
104 
105 // ------------- ----Copy constructor ---------------------------------
107  : FairField(),
108  fFileName(L.fFileName),
109  fPmin(L.fPmin), fPmax(L.fPmax), fPstep(L.fPstep),
110  fXmin(L.fXmin), fXmax(L.fXmax), fXstep(L.fXstep),
111  fYmin(L.fYmin), fYmax(L.fYmax), fYstep(L.fYstep),
112  fTmin(L.fTmin), fTmax(L.fTmax), fTstep(L.fTstep),
113  fFmin(L.fFmin), fFmax(L.fFmax), fFstep(L.fFstep),
114  fNp(L.fNp), fNx(L.fNx), fNy(L.fNy), fNt(L.fNt), fNf(L.fNf),
115  fBetaMean(L.fBetaMean), fBetaSig(L.fBetaSig), fBetaEff(L.fBetaEff)
116 {
117  fType = L.fType;
118 }
119 // ------------------------------------------------------------------------
120 
121 
122 
123 // ------------ Destructor --------------------------------------------
125  if ( fBetaMean ) delete fBetaMean;
126  if ( fBetaSig ) delete fBetaSig;
127  if ( fBetaEff ) delete fBetaEff;
128 }
129 // ------------------------------------------------------------------------
130 
131 
132 
133 // ----------- Intialisation ------------------------------------------
135  if (fFileName.EndsWith(".root")) ReadRootFile(fFileName, fName);
136  else if (fFileName.EndsWith(".dat")) ReadAsciiFile(fFileName);
137  else {
138  LOG(ERROR) << "-E- PndRichCalDb::Init: No proper file name defined! ("<<fFileName.Data()<<")";
139  Fatal("Init", "No proper file name");
140  }
141 }
142 // ------------------------------------------------------------------------
143 
144 
145 
146 // ----------- Get beta mean ---------------------------
148 
149  Double_t p = pnt.beta;
150  Double_t x = pnt.x;
151  Double_t y = pnt.y;
152  Double_t t = pnt.theta*180/3.1415927;
153  Double_t f = pnt.phi*180/3.1415927;
154  Int_t ip = 0;
155  Int_t ix = 0;
156  Int_t iy = 0;
157  Int_t it = 0;
158  Int_t iq = 0;
159  Double_t dp = 0.;
160  Double_t dx = 0.;
161  Double_t dy = 0.;
162  Double_t dt = 0.;
163  Double_t df = 0.;
164 
165  if (x>=0&&y<0) { f = 360 - f; y = -y; }
166  if (x<0&&y>=0) { f = 180 - f; x = -x; }
167  if (x<0&&y<0) { f = 180 + f; x = -x; y = -y; }
168  if (f>360) f -= 360;
169  if (f<0) f +=360;
170 
171  if ( IsInside(p, x, y, t, f, ip, ix, iy, it, iq, dp, dx, dy, dt, df) ) {
172 
173  Double_t pq[5][2];
174  pq[0][1] = dp; pq[0][0] = 1 - pq[0][1];
175  pq[1][1] = dx; pq[1][0] = 1 - pq[1][1];
176  pq[2][1] = dy; pq[2][0] = 1 - pq[2][1];
177  pq[3][1] = dt; pq[3][0] = 1 - pq[3][1];
178  pq[4][1] = df; pq[4][0] = 1 - pq[4][1];
179 
180  Double_t value = 0;
181  Int_t index;
182  for(Int_t bp=0; bp<2; bp++ ) {
183  for(Int_t bx=0; bx<2; bx++ ) {
184  for(Int_t by=0; by<2; by++ ) {
185  for(Int_t bt=0; bt<2; bt++ ) {
186  for(Int_t bf=0; bf<2; bf++ ) {
187  index =
188  (ip+bp)*fNx*fNy*fNt*fNf +
189  (ix+bx)*fNy*fNt*fNf +
190  (iy+by)*fNt*fNf +
191  (it+bt)*fNf +
192  (iq+bf);
193  value += fBetaMean->At(index)*pq[0][bp]*pq[1][bx]*pq[2][by]*pq[3][bt]*pq[4][bf];
194  }
195  }
196  }
197  }
198  }
199 
200  // Return interpolated field value
201  return value;
202 
203  }
204 
205  return 0.;
206 }
207 // ------------------------------------------------------------------------
208 
209 
210 
211 // ----------- Get beta sigma ---------------------------
213 
214  Double_t p = pnt.beta;
215  Double_t x = pnt.x;
216  Double_t y = pnt.y;
217  Double_t t = pnt.theta*180/3.1415927;
218  Double_t f = pnt.phi*180/3.1415927;
219  Int_t ip = 0;
220  Int_t ix = 0;
221  Int_t iy = 0;
222  Int_t it = 0;
223  Int_t iq = 0;
224  Double_t dp = 0.;
225  Double_t dx = 0.;
226  Double_t dy = 0.;
227  Double_t dt = 0.;
228  Double_t df = 0.;
229 
230  if (x>=0&&y<0) { f = 360 - f; y = -y; }
231  if (x<0&&y>=0) { f = 180 - f; x = -x; }
232  if (x<0&&y<0) { f = 180 + f; x = -x; y = -y; }
233  if (f>360) f -= 360;
234  if (f<0) f +=360;
235 
236  if ( IsInside(p, x, y, t, f, ip, ix, iy, it, iq, dp, dx, dy, dt, df) ) {
237 
238  Double_t pq[5][2];
239  pq[0][1] = dp; pq[0][0] = 1 - pq[0][1];
240  pq[1][1] = dx; pq[1][0] = 1 - pq[1][1];
241  pq[2][1] = dy; pq[2][0] = 1 - pq[2][1];
242  pq[3][1] = dt; pq[3][0] = 1 - pq[3][1];
243  pq[4][1] = df; pq[4][0] = 1 - pq[4][1];
244 
245  Double_t value = 0;
246  Int_t index;
247  for(Int_t bp=0; bp<2; bp++ ) {
248  for(Int_t bx=0; bx<2; bx++ ) {
249  for(Int_t by=0; by<2; by++ ) {
250  for(Int_t bt=0; bt<2; bt++ ) {
251  for(Int_t bf=0; bf<2; bf++ ) {
252  index =
253  (ip+bp)*fNx*fNy*fNt*fNf +
254  (ix+bx)*fNy*fNt*fNf +
255  (iy+by)*fNt*fNf +
256  (it+bt)*fNf +
257  (iq+bf);
258  value += fBetaSig->At(index)*pq[0][bp]*pq[1][bx]*pq[2][by]*pq[3][bt]*pq[4][bf];
259  }
260  }
261  }
262  }
263  }
264 
265  // Return interpolated field value
266  return value;
267 
268  }
269 
270  return 0.;
271 }
272 // ------------------------------------------------------------------------
273 
274 
275 // ----------- Get beta sigma ---------------------------
277 
278  Double_t p = pnt.beta;
279  Double_t x = pnt.x;
280  Double_t y = pnt.y;
281  Double_t t = pnt.theta*180/3.1415927;
282  Double_t f = pnt.phi*180/3.1415927;
283  Int_t ip = 0;
284  Int_t ix = 0;
285  Int_t iy = 0;
286  Int_t it = 0;
287  Int_t iq = 0;
288  Double_t dp = 0.;
289  Double_t dx = 0.;
290  Double_t dy = 0.;
291  Double_t dt = 0.;
292  Double_t df = 0.;
293 
294  if (x>=0&&y<0) { f = 360 - f; y = -y; }
295  if (x<0&&y>=0) { f = 180 - f; x = -x; }
296  if (x<0&&y<0) { f = 180 + f; x = -x; y = -y; }
297  if (f>360) f -= 360;
298  if (f<0) f +=360;
299 
300  if ( IsInside(p, x, y, t, f, ip, ix, iy, it, iq, dp, dx, dy, dt, df) ) {
301 
302  Double_t pq[5][2];
303  pq[0][1] = dp; pq[0][0] = 1 - pq[0][1];
304  pq[1][1] = dx; pq[1][0] = 1 - pq[1][1];
305  pq[2][1] = dy; pq[2][0] = 1 - pq[2][1];
306  pq[3][1] = dt; pq[3][0] = 1 - pq[3][1];
307  pq[4][1] = df; pq[4][0] = 1 - pq[4][1];
308 
309  Double_t value = 0;
310  Int_t index;
311  for(Int_t bp=0; bp<2; bp++ ) {
312  for(Int_t bx=0; bx<2; bx++ ) {
313  for(Int_t by=0; by<2; by++ ) {
314  for(Int_t bt=0; bt<2; bt++ ) {
315  for(Int_t bf=0; bf<2; bf++ ) {
316  index =
317  (ip+bp)*fNx*fNy*fNt*fNf +
318  (ix+bx)*fNy*fNt*fNf +
319  (iy+by)*fNt*fNf +
320  (it+bt)*fNf +
321  (iq+bf);
322  value += fBetaEff->At(index)*pq[0][bp]*pq[1][bx]*pq[2][by]*pq[3][bt]*pq[4][bf];
323  }
324  }
325  }
326  }
327  }
328 
329  // Return interpolated field value
330  return value;
331 
332  }
333 
334  return 0.;
335 }
336 // ------------------------------------------------------------------------
337 
338 
339 
340 
341 // ----------- Check whether a point is inside the map ----------------
343  Int_t& ip, Int_t& ix, Int_t& iy, Int_t& it, Int_t& iq,
344  Double_t& dp, Double_t& dx, Double_t& dy, Double_t& dt, Double_t& df) {
345 
346  // --- Check for being outside the map range
347  p = p < fPmax ? p : fPmax-0.00001;
348  t = t < fTmax ? t : fTmax-0.01;
349  if ( ! ( p >= fPmin && p < fPmax &&
350  x >= fXmin && x < fXmax && y >= fYmin && y < fYmax &&
351  t >= fTmin && t < fTmax && f >= fFmin && f < fFmax) ) {
352  ip = ix = iy = it = iq = 0;
353  dp = dx = dy = dt = df = 0.;
354  return kFALSE;
355  }
356 
357  dp = ( p - fPmin ) / fPstep;
358  dx = ( x - fXmin ) / fXstep;
359  dy = ( y - fYmin ) / fYstep;
360  dt = ( t - fTmin ) / fTstep;
361  df = ( f - fFmin ) / fFstep;
362 
363  // --- Determine grid cell
364  ip = Int_t( dp );
365  ix = Int_t( dx );
366  iy = Int_t( dy );
367  it = Int_t( dt );
368  iq = Int_t( df );
369 
370 
371  // Relative distance from grid point (in units of cell size)
372  dp -= Double_t(ip);
373  dx -= Double_t(ix);
374  dy -= Double_t(iy);
375  dt -= Double_t(it);
376  df -= Double_t(iq);
377 
378  return kTRUE;
379 
380 }
381 // ------------------------------------------------------------------------
382 
383 
384 
385 // ---------- Write the map to an ASCII file --------------------------
386 void PndRichCalDb::WriteAsciiFile(const char* fileName) {
387 
388  // Open file
389  LOG(INFO) << "PndRichCalDb: Writing field map to ASCII file "<<fileName;
390  ofstream mapFile(fileName);
391  if ( ! mapFile.is_open() ) {
392  LOG(ERROR) << "PndRichCalDb:ReadAsciiFile: Could not open file! ";
393  return;
394  }
395 
396  // Write field map grid parameters
397  mapFile.precision(4);
398  mapFile << showpoint;
399 
400  mapFile << fPmin << " " << fPmax << " " << fNp << endl;
401  mapFile << fXmin << " " << fXmax << " " << fNx << endl;
402  mapFile << fYmin << " " << fYmax << " " << fNy << endl;
403  mapFile << fTmin << " " << fTmax << " " << fNt << endl;
404  mapFile << fFmin << " " << fFmax << " " << fNf << endl;
405 
406  // Write field values
407  cout << right;
408  //Int_t nTot = fNp * fNx * fNy * fNt * fNf; //[R.K. 01/2017] unused variable?
409  cout << "-I- PndRichCalDb: " << fNp*fNx*fNy*fNt*fNf << " entries to write... "
410  << setw(3) << 0 << " % ";
411  Int_t index=0;
412 // div_t modul;
413 // Int_t iDiv = TMath::Nint(nTot/100.);
414  for(Int_t ip=0; ip<fNp; ip++) {
415  for(Int_t ix=0; ix<fNx; ix++) {
416  for(Int_t iy=0; iy<fNy; iy++) {
417  for(Int_t it=0; it<fNt; it++) {
418  for(Int_t iq=0; iq<fNf; iq++) {
419  index =
420  ip*fNx*fNy*fNt*fNf +
421  ix*fNy*fNt*fNf +
422  iy*fNt*fNf +
423  it*fNf +
424  iq;
425 /* if (iDiv!=0)
426  {
427  modul = div(index,iDiv);
428  if ( modul.rem == 0 ) {
429  Double_t perc = TMath::Nint(100.*index/nTot);
430  cout << "\b\b\b\b\b\b" << setw(3) << perc << " % " << flush;
431  }
432  }*/
433  mapFile << fBetaMean->At(index) << " " << fBetaSig->At(index)
434  << " " << fBetaEff->At(index) << endl;
435  } // f-Loop
436  } // t-Loop
437  } // y-Loop
438  } // x-Loop
439  } // p-Loop
440  mapFile.close();
441 
442 }
443 // ------------------------------------------------------------------------
444 
445 
446 
447 // ------- Write field map to a ROOT file -----------------------------
448 void PndRichCalDb::WriteRootFile(const char* fileName,
449  const char* mapName) {
450 
451  PndRichCalDbData* data = new PndRichCalDbData(mapName, *this);
452  TFile* oldFile = gFile;
453  TFile* file = new TFile(fileName, "RECREATE");
454  data->Write();
455  file->Close();
456  if(oldFile) oldFile->cd();
457 
458 }
459 // ------------------------------------------------------------------------
460 
461 
462 
463 // --------- Screen output --------------------------------------------
465  TString type = "Map";
466  if ( fType == 2 ) type = "Soleniod Map ";
467  if ( fType == 3 ) type = "Dipole Map ";
468  if ( fType == 4 ) type = "Trans Map ";
469  cout << "======================================================" << endl;
470  cout.precision(4);
471  cout << showpoint;
472  cout << "---- " << fTitle << " : " << fName << endl;
473  cout << "----" << endl;
474  cout << "---- Field type : " << type << endl;
475  cout << "----" << endl;
476  cout << "---- Field map grid : " << endl;
477  cout << "---- p = " << setw(4) << fPmin << " to " << setw(4) << fPmax
478  << " , " << fNp << " grid points, dp = " << fPstep << " " << endl;
479  cout << "---- x = " << setw(4) << fXmin << " to " << setw(4) << fXmax
480  << " cm, " << fNx << " grid points, dx = " << fXstep << " cm" << endl;
481  cout << "---- y = " << setw(4) << fYmin << " to " << setw(4) << fYmax
482  << " cm, " << fNy << " grid points, dy = " << fYstep << " cm" << endl;
483  cout << "---- t = " << setw(4) << fTmin << " to " << setw(4) << fTmax
484  << " , " << fNt << " grid points, dt = " << fTstep << " " << endl;
485  cout << "---- f = " << setw(4) << fFmin << " to " << setw(4) << fFmax
486  << " , " << fNf << " grid points, df = " << fFstep << " " << endl;
487  cout << endl;
488 /* cout << "---- Field scaling factor: " << fScale << endl;
489  Double_t bx = GetBx(0.,0.,0.);
490  Double_t by = GetBy(0.,0.,0.);
491  Double_t bz = GetBz(0.,0.,0.);
492  cout << "----" << endl;
493  cout << "---- Field at origin is ( " << setw(6) << bx << ", " << setw(6)
494  << by << ", " << setw(6) << bz << ") kG" << endl;*/
495  cout << "======================================================" << endl;
496 }
497 // ------------------------------------------------------------------------
498 
499 
500 
501 // --------- Reset parameters and data (private) ----------------------
503  fPmin = fXmin = fYmin = fTmin = fFmin = 0;
504  fPmax = fXmax = fYmax = fTmax = fFmax = 0;
505  fPstep = fXstep = fYstep = fTstep = fFstep = 0;
506  fNp = fNx = fNy = fNt = fNf = 0;
507  if ( fBetaMean ) { delete fBetaMean; fBetaMean = NULL; }
508  if ( fBetaSig ) { delete fBetaSig; fBetaSig = NULL; }
509  if ( fBetaEff ) { delete fBetaEff; fBetaEff = NULL; }
510 }
511 // ------------------------------------------------------------------------
512 
513 
514 
515 // ----- Read field map from ASCII file (private) ---------------------
516 void PndRichCalDb::ReadAsciiFile(const char* fileName) {
517 
518  // Open file
519  cout << "-I- PndRichCalDb: Reading field map from ASCII file "
520  << fileName << endl;
521  ifstream mapFile(fileName);
522  if ( ! mapFile.is_open() ) {
523  cerr << "-E- PndRichCalDb:ReadAsciiFile: Could not open file! " << endl;
524  Fatal("ReadAsciiFile","Could not open file");
525  }
526 
527  // Read grid parameters
528 
529  mapFile >>fPmin >> fPmax >> fNp;
530  mapFile >>fXmin >> fXmax >> fNx;
531  mapFile >>fYmin >> fYmax >> fNy;
532  mapFile >>fTmin >> fTmax >> fNt;
533  mapFile >>fFmin >> fFmax >> fNf;
534  fPstep = ( fPmax - fPmin ) / Double_t( fNp - 1 );
535  fXstep = ( fXmax - fXmin ) / Double_t( fNx - 1 );
536  fYstep = ( fYmax - fYmin ) / Double_t( fNy - 1 );
537  fTstep = ( fTmax - fTmin ) / Double_t( fNt - 1 );
538  fFstep = ( fFmax - fFmin ) / Double_t( fNf - 1 );
539 
540  // Create field arrays
541  fBetaMean = new TArrayF(fNp * fNx * fNy * fNt * fNf);
542  fBetaSig = new TArrayF(fNp * fNx * fNy * fNt * fNf);
543  fBetaEff = new TArrayF(fNp * fNx * fNy * fNt * fNf);
544 
545  // Read the field values
546  cout << right;
547  Int_t nTot = fNp * fNx * fNy * fNt * fNf;
548  cout << "-I- PndRichCalDb: " << nTot << " entries to read... "
549  << setw(3) << 0 << " % ";
550  Int_t index = 0;
551 // div_t modul;
552 // Int_t iDiv = TMath::Nint(nTot/100.);
553  Double_t m, s, e;
554  for (Int_t ip=0; ip<fNp; ip++) {
555  for (Int_t ix=0; ix<fNx; ix++) {
556  for (Int_t iy = 0; iy<fNy; iy++) {
557  for (Int_t it = 0; it<fNt; it++) {
558  for (Int_t iq = 0; iq<fNf; iq++) {
559  if (! mapFile.good()) cerr << "-E- PndRichCalDb::ReadAsciiFile: "
560  << "I/O Error at " << ip << " " << ix << " "
561  << iy << " " << it << " " << iq << endl;
562  index =
563  ip*fNx*fNy*fNt*fNf +
564  ix*fNy*fNt*fNf +
565  iy*fNt*fNf +
566  it*fNf +
567  iq;
568 /* if (iDiv!=0)
569  {
570  modul = div(index,iDiv);
571  if ( modul.rem == 0 ) {
572  Double_t perc = TMath::Nint(100.*index/nTot);
573  cout << "\b\b\b\b\b\b" << setw(3) << perc << " % " << flush;
574  }
575  }*/
576  mapFile >> m >> s >> e;
577  if (s<0||e<0||e>1)
578  cout << index << " " << m << " " << s << " " << e << endl;
579  fBetaMean->AddAt(m, index);
580  fBetaSig->AddAt(s, index);
581  fBetaEff->AddAt(e, index);
582  if ( mapFile.eof() ) {
583  cerr << endl << "-E- PndRichCalDb::ReadAsciiFile: EOF"
584  << " reached at " << ip << " " << ix << " " << iy << " " << it << " " << iq << endl;
585  mapFile.close();
586  break;
587  }
588  } // f-Loop
589  } // t-Loop
590  } // y-Loop
591  } // x-Loop
592  } // p-Loop
593 
594  cout << " " << index+1 << " read" << endl;
595 
596  mapFile.close();
597 
598 }
599 // ------------------------------------------------------------------------
600 
601 
602 
603 // ------------- Read field map from ROOT file (private) ---------------
604 void PndRichCalDb::ReadRootFile(const char* fileName,
605  const char* mapName) {
606 
607  // Store gFile pointer
608  TFile* oldFile = gFile;
609 
610  // Open root file
611  LOG(INFO) << "PndRichCalDb: Reading field map from ROOT file "<<fileName;
612  TFile* file = new TFile(fileName, "READ");
613  if (file->IsZombie()) {
614  LOG(ERROR) << "-E- PndRichCalDb::ReadRootfile: Cannot read from file! (" << fileName <<")";
615  Fatal("ReadRootFile","Cannot read from file");
616  }
617 
618  // Get the field data object
619  PndRichCalDbData* data = NULL;
620  file->GetObject(mapName, data);
621  if ( ! data ) {
622  LOG(ERROR) << "PndRichCalDb::ReadRootFile: data object " << fileName << " not found in file! ";
623  exit(-1);
624  }
625 
626  // Get the field parameters
627  SetCalDb(data);
628 
629  // Close the root file and delete the data object
630  file->Close();
631  delete data;
632  delete file;
633  if ( oldFile ) oldFile->cd();
634 
635 }
636 // ------------------------------------------------------------------------
637 
638 
639 
640 // ------------ Set field parameters and data (private) ----------------
642 
643  // Check compatibility
644  std::cout << "fType = " << fType << " " << data->GetType() << std::endl;
645  fType = data->GetType();
646  if ( data->GetType() != fType ) {
647  LOG(ERROR) << "PndRichCalDb::SetField: Incompatible map types Field map is of type "<<fType<<" \n but map on file is of type " << data->GetType();
648  Fatal("SetField","Incompatible map types");
649  }
650 
651 
652  fPmin = data->GetPmin();
653  fXmin = data->GetXmin();
654  fYmin = data->GetYmin();
655  fTmin = data->GetTmin();
656  fFmin = data->GetFmin();
657  fPmax = data->GetPmax();
658  fXmax = data->GetXmax();
659  fYmax = data->GetYmax();
660  fTmax = data->GetTmax();
661  fFmax = data->GetFmax();
662  fNp = data->GetNp();
663  fNx = data->GetNx();
664  fNy = data->GetNy();
665  fNt = data->GetNt();
666  fNf = data->GetNf();
667  fPstep = ( fPmax - fPmin ) / Double_t( fNp - 1 );
668  fXstep = ( fXmax - fXmin ) / Double_t( fNx - 1 );
669  fYstep = ( fYmax - fYmin ) / Double_t( fNy - 1 );
670  fTstep = ( fTmax - fTmin ) / Double_t( fNt - 1 );
671  fFstep = ( fFmax - fFmin ) / Double_t( fNf - 1 );
672  if ( fBetaMean ) delete fBetaMean;
673  if ( fBetaSig ) delete fBetaSig;
674  if ( fBetaEff ) delete fBetaEff;
675  fBetaMean = new TArrayF(*(data->GetBetaMean()));
676  fBetaSig = new TArrayF(*(data->GetBetaSig()));
677  fBetaEff = new TArrayF(*(data->GetBetaEff()));
678 
679 }
680 // ------------------------------------------------------------------------
681 
682 
683 
684 
Double_t GetTmax() const
Int_t GetNt() const
virtual void Print()
Double_t p
Definition: anasim.C:58
Double_t x
Definition: PndRichCalDb.h:25
double dy
TArrayF * fBetaSig
Definition: PndRichCalDb.h:165
cout<< "-----------------------------------------------> Quarter VOLUME<<endl;name="QuarterShape";QuarterShape=newTGeoArb8(name,dz,vertQuar);name="Quarter4Vol";TStringmedium="air";QuarterVol=newTGeoVolumeAssembly(name);name="SubunitShape";SubunitShape=newTGeoArb8(name,dz,vertSub);TStringmedium="air";name="SubunitVol";name1="SubunitVol1";name2="SubunitVol2";name3="SubunitVol3";name4="SubunitVol4";name5="SubunitVol5";name6="SubunitVol6";name7="SubunitVol7";name8="SubunitVol8";name9="SubunitVol9";SubunitVol=newTGeoVolumeAssembly(name);SubunitVol1=newTGeoVolumeAssembly(name1);SubunitVol2=newTGeoVolumeAssembly(name2);SubunitVol3=newTGeoVolumeAssembly(name3);SubunitVol4=newTGeoVolumeAssembly(name4);SubunitVol5=newTGeoVolumeAssembly(name5);SubunitVol6=newTGeoVolumeAssembly(name6);SubunitVol7=newTGeoVolumeAssembly(name7);SubunitVol8=newTGeoVolumeAssembly(name8);SubunitVol9=newTGeoVolumeAssembly(name9);name="BoxShape";BoxShape=newTGeoArb8(name,dz,vertBox);TStringmedium="air";name="BoxVol";BoxVol=newTGeoVolumeAssembly(name);name1="BoxVol1";name2="BoxVol2";name3="BoxVol3";name4="BoxVol4";BoxVol1=newTGeoVolumeAssembly(name1);BoxVol2=newTGeoVolumeAssembly(name2);BoxVol3=newTGeoVolumeAssembly(name3);BoxVol4=newTGeoVolumeAssembly(name4);for(Int_tb=0;b<kNumOfBoxes;b++){cout<<""<<endl;cout<<"---------------->BOXnumber:"<<b<<endl;if(b==0){trBox=newTGeoTranslation(tr+kSpaceInSub+kAlveoleThickness+0.5*kSpaceInBox,tr+kSpaceInSub+kAlveoleThickness+0.5*kSpaceInBox,0.);rotBox=newTGeoRotation();rotBox.RotateX(0.465518);rotBox.RotateY(-0.465518);}if(b==1){trBox=newTGeoTranslation(-tr-kSpaceInSub-kAlveoleThickness-0.5*kSpaceInBox,-tr-kSpaceInSub-kAlveoleThickness-0.5*kSpaceInBox,0.);rotBox=newTGeoRotation();rotBox.RotateX(-0.465518);rotBox.RotateY(0.465518);}if(b==2){trBox=newTGeoTranslation(tr+kSpaceInSub+kAlveoleThickness+0.5*kSpaceInBox,-tr-kSpaceInSub-kAlveoleThickness-0.5*kSpaceInBox,0.);rotBox=newTGeoRotation();rotBox.RotateX(-0.465518);rotBox.RotateY(-0.465518);}if(b==3){trBox=newTGeoTranslation(-tr-kSpaceInSub-kAlveoleThickness-0.5*kSpaceInBox,tr+kSpaceInSub+kAlveoleThickness+0.5*kSpaceInBox,0.);rotBox=newTGeoRotation();rotBox.RotateX(0.465518);rotBox.RotateY(0.465518);}TGeoCombiTrans*trrotBox=newTGeoCombiTrans(trBox,rotBox);name="BoxVol";name+=b;trrotBox->SetName(name);trrotBox->RegisterYourself();SubunitVol->AddNode(BoxVol,b,trrotBox);if(b==1){name+=b;trrotBox->SetName(name);trrotBox->RegisterYourself();SubunitVol1->AddNode(BoxVol,b,trrotBox);}if(b==2){name+=b;trrotBox->SetName(name);trrotBox->RegisterYourself();SubunitVol2->AddNode(BoxVol1,b,trrotBox);}if(b==0){name+=b;trrotBox-> SetName(name)
Double_t fPstep
Definition: PndRichCalDb.h:152
Double_t GetFmax() const
__m128 m
Definition: P4_F32vec4.h:28
TArrayF * fBetaMean
Definition: PndRichCalDb.h:164
TFile * file
Double_t GetYmax() const
TString fFileName
Definition: PndRichCalDb.h:148
Double_t fTmax
Definition: PndRichCalDb.h:155
exit(0)
Int_t GetType() const
Double_t fYstep
Definition: PndRichCalDb.h:154
TLorentzVector s
Definition: Pnd2DStar.C:50
Double_t fFstep
Definition: PndRichCalDb.h:156
TClonesArray * pnt
TArrayF * GetBetaSig() const
Int_t GetType() const
void ReadRootFile(const char *fileName, const char *mapName)
TArrayF * fBetaEff
Definition: PndRichCalDb.h:166
Double_t fXstep
Definition: PndRichCalDb.h:153
TArrayF * GetBetaSig() const
Definition: PndRichCalDb.h:112
Int_t GetNf() const
void ReadAsciiFile(const char *fileName)
Double_t fTstep
Definition: PndRichCalDb.h:155
TArrayF * GetBetaMean() const
Definition: PndRichCalDb.h:111
Int_t GetNp() const
TArrayF * GetBetaEff() const
Definition: PndRichCalDb.h:113
Double_t GetXmin() const
Double_t fTmin
Definition: PndRichCalDb.h:155
Double_t fYmax
Definition: PndRichCalDb.h:154
Int_t GetNx() const
void WriteAsciiFile(const char *fileName)
Double_t fPmin
Definition: PndRichCalDb.h:152
Double_t GetTmin() const
Double_t GetYmin() const
TArrayF * GetBetaMean() const
Double_t
virtual Bool_t IsInside(Double_t p, Double_t x, Double_t y, Double_t t, Double_t f, Int_t &ip, Int_t &ix, Int_t &iy, Int_t &it, Int_t &iq, Double_t &dp, Double_t &dx, Double_t &dy, Double_t &dt, Double_t &df)
void WriteRootFile(const char *fileName, const char *mapName)
void MapName(TString &name)
TArrayF * GetBetaEff() const
TFile * f
Definition: bump_analys.C:12
Double_t fPmax
Definition: PndRichCalDb.h:152
Double_t fFmax
Definition: PndRichCalDb.h:156
double dx
void SetCalDb(const PndRichCalDbData *data)
Double_t fYmin
Definition: PndRichCalDb.h:154
Double_t x
ClassImp(PndAnaContFact)
Double_t theta
Definition: PndRichCalDb.h:27
Double_t fXmax
Definition: PndRichCalDb.h:153
TTree * t
Definition: bump_analys.C:13
virtual void Init()
Double_t phi
Definition: PndRichCalDb.h:28
Double_t fXmin
Definition: PndRichCalDb.h:153
Double_t y
TF1 * bp
Definition: hist-t7.C:69
virtual ~PndRichCalDb()
Double_t GetFmin() const
Double_t beta
Definition: PndRichCalDb.h:24
Double_t GetPmax() const
Int_t GetNy() const
Double_t GetPmin() const
Double_t y
Definition: PndRichCalDb.h:26
Double_t fFmin
Definition: PndRichCalDb.h:156
Double_t GetXmax() const