/[PAMELA software]/DarthVader/CalorimeterLevel2/src/CaloLevel0.cpp
ViewVC logotype

Contents of /DarthVader/CalorimeterLevel2/src/CaloLevel0.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.1 - (show annotations) (download)
Sun Mar 18 20:22:54 2007 UTC (17 years, 8 months ago) by mocchiut
Branch: MAIN
Renamed CaloProcessing into CaloLevel0, added method to calorimeter level0 retrieving, fixed OrbitalInfo class doxygen bug

1 /**
2 * \file src/CaloLevel0.cpp
3 * \author Emiliano Mocchiutti
4 **/
5 //
6 // C/C++ headers
7 //
8 #include <sstream>
9 #include <fstream>
10 //
11 // ROOT headers
12 //
13 #include <TTree.h>
14 #include <TBranch.h>
15 #include <TFile.h>
16 #include <TObject.h>
17 //
18 // YODA headers
19 //
20 #include <PamelaRun.h>
21 #include <physics/calorimeter/CalorimeterEvent.h>
22 #include <CalibCalPedEvent.h>
23 //
24 //
25 //
26 #include <GLTables.h>
27 //
28 // this package headers
29 //
30 #include <delay.h>
31 #include <CaloLevel0.h>
32 //
33 //
34 // Declaration of the core fortran routines
35 //
36 #define calol2cm calol2cm_
37 extern "C" int calol2cm();
38 #define calol2tr calol2tr_
39 extern "C" int calol2tr();
40 //
41 using namespace std;
42 //
43 //
44 // Public methods
45 //
46
47 CaloLevel0::~CaloLevel0(){
48 if ( de ) delete de;
49 delete this;
50 }
51
52 CaloLevel0::CaloLevel0(){
53 //
54 extern struct FlCaLevel1 clevel1_;
55 extern struct FlCaLevel2 clevel2_;
56 clevel1 = &clevel1_;
57 clevel2 = &clevel2_;
58 //
59 trkseqno = 0;
60 ClearStructs();
61 //
62 memset(dexy, 0, 2*22*96*sizeof(Float_t));
63 memset(dexyc, 0, 2*22*96*sizeof(Float_t));
64 memset(mip, 0, 2*22*96*sizeof(Float_t));
65 memset(base, 0, 2*22*6*sizeof(Float_t));
66 memset(sbase, 0, 2*22*6*sizeof(Float_t));
67 calopar1 = true;
68 calopar2 = true;
69 calopar3 = true;
70 crosst = true;
71 ftcalopar1 = 0;
72 ttcalopar1 = 0;
73 ftcalopar2 = 0;
74 ttcalopar2 = 0;
75 ftcalopar3 = 0;
76 ttcalopar3 = 0;
77 }
78
79 void CaloLevel0::SetCrossTalk(Bool_t ct){
80 crosst = ct;
81 };
82
83 void CaloLevel0::SetVerbose(Bool_t ct){
84 verbose = ct;
85 };
86
87 /**
88 * Initialize CaloLevel0 object
89 **/
90 void CaloLevel0::ProcessingInit(TSQLServer *dbc, UInt_t hs, Int_t &sgnl, TTree *l0tree, Bool_t isdeb, Bool_t isverb){
91 //
92 debug = isdeb;
93 verbose = isverb;
94 //
95 l0tr=(TTree*)l0tree;
96 de = new pamela::calorimeter::CalorimeterEvent();
97 l0calo = (TBranch*)l0tr->GetBranch("Calorimeter");
98 l0tr->SetBranchAddress("Calorimeter", &de);
99 //
100 trkseqno = 0;
101 ClearStructs();
102 //
103 GL_CALO_CALIB *glcalo = new GL_CALO_CALIB();
104 //
105 sgnl = 0;
106 UInt_t uptime = 0;
107 //
108 for (Int_t s = 0; s < 4; s++){
109 idcalib[s] = 0;
110 fromtime[s] = 0;
111 totime[s] = 0;
112 calibno[s] = 0;
113 ClearCalibVals(s);
114 //
115 sgnl = glcalo->Query_GL_CALO_CALIB(hs,uptime,s,dbc);
116 if ( sgnl < 0 ){
117 if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
118 return;
119 };
120 //
121 idcalib[s] = glcalo->ID_ROOT_L0;
122 fromtime[s] = glcalo->FROM_TIME;
123 if ( glcalo->TO_TIME < hs ){ // calibration is corrupted and we are using the one that preceed the good one
124 totime[s] = uptime;
125 } else {
126 totime[s] = glcalo->TO_TIME;
127 };
128 calibno[s] = glcalo->EV_ROOT;
129 //
130 if ( totime[s] == 0 ){
131 if ( verbose ) printf(" CALORIMETER - WARNING: data with no associated calibration\n");
132 ClearCalibVals(s);
133 sgnl = 100;
134 };
135 };
136 //
137 // determine path and name and entry of the calibration file
138 //
139 GL_ROOT *glroot = new GL_ROOT();
140 if ( verbose ) printf("\n");
141 for (Int_t s = 0; s < 4; s++){
142 if ( verbose ) printf(" ** SECTION %i **\n",s);
143 if ( totime[s] > 0 ){
144 //
145 sgnl = glroot->Query_GL_ROOT(idcalib[s],dbc);
146 if ( sgnl < 0 ){
147 if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
148 return;
149 };
150 //
151 stringstream name;
152 name.str("");
153 name << glroot->PATH.Data() << "/";
154 name << glroot->NAME.Data();
155 //
156 fcalname[s] = (TString)name.str().c_str();
157 if ( verbose ) printf(" - runheader at time %u. From time %u to time %u \n use file %s \n calibration at entry %i \n\n",hs,fromtime[s],totime[s],fcalname[s].Data(),calibno[s]);
158 } else {
159 if ( verbose ) printf(" - runheader at time %u. NO CALIBRATION INCLUDE THE RUNHEADER! ",hs);
160 };
161 sgnl = LoadCalib(s);
162 if ( sgnl ) break;
163 };
164 //
165 delete glcalo;
166 delete glroot;
167 //
168 return;
169 //
170 }
171
172 Int_t CaloLevel0::ChkCalib(TSQLServer *dbc, UInt_t atime){
173 Int_t sgnl = 0;
174 for ( Int_t s = 0; s < 4; s++){
175 if ( atime > totime[s] ){
176 if ( !dbc->IsConnected() ) throw -116;
177 sgnl = Update(dbc,atime,s);
178 if ( sgnl < 0 ) return(sgnl);
179 };
180 };
181 return(sgnl);
182 }
183
184 Int_t CaloLevel0::ChkParam(TSQLServer *dbc, UInt_t runheader){
185 stringstream calfile;
186 stringstream bmfile;
187 stringstream aligfile;
188 Int_t error = 0;
189 FILE *f = 0;
190 ifstream badfile;
191 GL_PARAM *glparam = new GL_PARAM();
192 //
193 if ( calopar1 || ( ttcalopar1 != 0 && ttcalopar1 < runheader ) ){
194 calopar1 = false;
195 //
196 // determine where I can find calorimeter ADC to MIP conversion file
197 //
198 if ( verbose ) printf(" Querying DB for calorimeter parameters files...\n");
199 //
200 error = 0;
201 error = glparam->Query_GL_PARAM(runheader,101,dbc);
202 if ( error < 0 ) return(error);
203 //
204 calfile.str("");
205 calfile << glparam->PATH.Data() << "/";
206 calfile << glparam->NAME.Data();
207 ftcalopar1 = glparam->FROM_TIME;
208 ttcalopar1 = glparam->TO_TIME;
209 //
210 if ( verbose ) printf("\n Using ADC to MIP conversion file: \n %s \n",calfile.str().c_str());
211 f = fopen(calfile.str().c_str(),"rb");
212 if ( !f ){
213 if ( verbose ) printf(" CALORIMETER - ERROR: no ADC to MIP file!\n");
214 return(-105);
215 };
216 //
217 for (Int_t m = 0; m < 2 ; m++ ){
218 for (Int_t k = 0; k < 22; k++ ){
219 for (Int_t l = 0; l < 96; l++ ){
220 fread(&mip[m][k][l],sizeof(mip[m][k][l]),1,f);
221 };
222 };
223 };
224 fclose(f);
225 };
226 //
227 if ( calopar2 || ( ttcalopar2 != 0 && ttcalopar2 < runheader ) ){
228 calopar2 = false;
229 //
230 // determine where I can find calorimeter alignment file
231 //
232 //
233 error = 0;
234 error = glparam->Query_GL_PARAM(runheader,102,dbc);
235 if ( error < 0 ) return(error);
236 //
237 aligfile.str("");
238 aligfile << glparam->PATH.Data() << "/";
239 aligfile << glparam->NAME.Data();
240 ftcalopar2 = glparam->FROM_TIME;
241 ttcalopar2 = glparam->TO_TIME;
242 //
243 if ( verbose ) printf("\n Using alignment file: \n %s \n",aligfile.str().c_str());
244 f = fopen(aligfile.str().c_str(),"rb");
245 if ( !f ){
246 if ( verbose ) printf(" CALORIMETER - ERROR: no alignement file!\n");
247 return(-106);
248 };
249 //
250 fread(&clevel1->xalig,sizeof(clevel1->xalig),1,f);
251 if ( debug ) printf(" xalig = %f \n",clevel1->xalig);
252 fread(&clevel1->yalig,sizeof(clevel1->yalig),1,f);
253 if ( debug ) printf(" yalig = %f \n",clevel1->yalig);
254 fread(&clevel1->zalig,sizeof(clevel1->zalig),1,f);
255 if ( debug ) printf(" zalig = %f \n",clevel1->zalig);
256 fread(&clevel1->emin,sizeof(clevel1->emin),1,f);
257 if ( debug ) printf(" signal threshold = %f \n",clevel1->emin);
258 //
259 fclose(f);
260 };
261 //
262 // Load offline bad strip mask
263 //
264 if ( calopar3 || ( ttcalopar3 != 0 && ttcalopar3 < runheader ) ){
265 calopar3 = false;
266 //
267 // determine where I can find calorimeter alignment file
268 //
269 //
270 error = 0;
271 error = glparam->Query_GL_PARAM(runheader,103,dbc);
272 if ( error < 0 ) return(error);
273 //
274 bmfile.str("");
275 bmfile << glparam->PATH.Data() << "/";
276 bmfile << glparam->NAME.Data();
277 ftcalopar3 = glparam->FROM_TIME;
278 ttcalopar3 = glparam->TO_TIME;
279 //
280 if ( verbose ) printf("\n Using bad strip offline mask file: \n %s \n\n",bmfile.str().c_str());
281 badfile.open(bmfile.str().c_str());
282 if ( !badfile ){
283 if ( verbose ) printf(" CALORIMETER - ERROR: no bad strip offline mask file!\n");
284 return(-115);
285 };
286 //
287 Bool_t isdone = false;
288 Int_t bad = 0;
289 Int_t view = 1;
290 Int_t strip = 0;
291 Int_t plane = 21;
292 while ( !isdone ) {
293 badfile >> bad;
294 obadmask[view][plane][strip] = bad;
295 if ( debug && bad ) printf(" SETTING view %i plane %i strip %i BAD = %i \n",view,plane,strip,bad);
296 strip++;
297 if ( strip > 95 ){
298 strip = 0;
299 plane--;
300 if ( plane < 0 ){
301 plane = 21;
302 view--;
303 };
304 if ( view < 0 ) isdone = true;
305 };
306 };
307 //
308 badfile.close();
309 };
310 //
311 delete glparam;
312 //
313 return(0);
314 }
315
316
317
318 void CaloLevel0::FindBaseRaw(Int_t l, Int_t m, Int_t pre){
319 Float_t minstrip = 100000.;
320 Float_t rms = 0.;
321 base[l][m][pre] = 0.;
322 for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
323 if ( calgood[l][m][e] == 0. && obadmask[l][m][e] == 0 && dexy[l][m][e]-calped[l][m][e] < minstrip && dexy[l][m][e] > 0.) {
324 minstrip = dexy[l][m][e]-calped[l][m][e];
325 rms = calthr[l][m][pre];
326 };
327 };
328 if ( minstrip != 100000. ) {
329 Float_t strip6s = 0.;
330 for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
331 if ( (dexy[l][m][e]-calped[l][m][e]) >= minstrip && (dexy[l][m][e]-calped[l][m][e]) <= (minstrip+rms) ) {
332 strip6s += 1.;
333 base[l][m][pre] += (dexy[l][m][e] - calped[l][m][e]);
334 };
335 //
336 // compression
337 //
338 if ( abs((int)(dexy[l][m][e]-calped[l][m][e])) <= (minstrip+rms) ) {
339 dexyc[l][m][e] = 0.;
340 } else {
341 dexyc[l][m][e] = dexy[l][m][e];
342 };
343 };
344 if ( strip6s >= 9. ){
345 Double_t arro = base[l][m][pre]/strip6s;
346 Float_t deci = 1000.*((float)arro - float(int(arro)));
347 if ( deci < 500. ) {
348 arro = double(int(arro));
349 } else {
350 arro = 1. + double(int(arro));
351 };
352 base[l][m][pre] = arro;
353 } else {
354 base[l][m][pre] = 31000.;
355 for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
356 dexyc[l][m][e] = dexy[l][m][e];
357 };
358 };
359 } else {
360 base[l][m][pre] = 31000.;
361 };
362 }
363
364 Int_t CaloLevel0::Calibrate(Int_t ei){
365 //
366 // get entry ei
367 //
368 l0calo->GetEntry(ei);
369 //
370 // if it was not a selftrigger event, could it ever been a selftrigger event? if so trigty = 3.
371 //
372 Int_t val = 0;
373 Int_t del = 1100;
374 if ( clevel2->trigty != 2. ){
375 Bool_t ck = false;
376 for (Int_t sec = 0; sec < 4; sec++){
377 val = (Int_t)de->calselftrig[sec][6];
378 del = delay(val);
379 if ( del < 1100 ){
380 clevel2->wartrig = 0.;
381 clevel2->trigty = 3.;
382 ck = true;
383 break;
384 };
385 };
386 if ( !ck ) clevel2->wartrig = 100.;
387 } else {
388 Bool_t ck = false;
389 for (Int_t sec = 0; sec < 4; sec++){
390 val = (Int_t)de->calselftrig[sec][6];
391 del = delay(val);
392 if ( del < 1100 ){
393 clevel2->wartrig = 0.;
394 ck = true;
395 };
396 };
397 if ( !ck ) clevel2->wartrig = 100.;
398 };
399 //
400 Int_t se = 5;
401 Int_t done = 0;
402 Int_t pre = -1;
403 Bool_t isCOMP = false;
404 Bool_t isFULL = false;
405 Bool_t isRAW = false;
406 Float_t ener;
407 Int_t doneb = 0;
408 Int_t donec = 0;
409 Int_t ck = 0;
410 Int_t ipre = 0;
411 Int_t ip[3] = {0};
412 Float_t base0, base1, base2;
413 base0 = 0.;
414 base1 = 0.;
415 base2 = 0.;
416 Float_t qpre[6] = {0.,0.,0.,0.,0.,0.};
417 Float_t ene[96];
418 Int_t chdone[4] = {0,0,0,0};
419 Int_t pe = 0;
420 //
421 Float_t ener0 = 0.;
422 Float_t cbase0 = 0.;
423 Bool_t pproblem = false;
424 //
425 Float_t tim = 0.;
426 Int_t plo = 0;
427 Int_t fbi = 0;
428 Int_t cle = 0;
429 //
430 // run over views and planes
431 //
432 for (Int_t l = 0; l < 2; l++){
433 for (Int_t m = 0; m < 22; m++){
434 //
435 // determine the section number
436 //
437 se = 5;
438 if (l == 0 && m%2 == 0) se = 3;
439 if (l == 0 && m%2 != 0) se = 2;
440 if (l == 1 && m%2 == 0) se = 1;
441 if (l == 1 && m%2 != 0) se = 0;
442 //
443 // determine what kind of event we are going to analyze
444 //
445 isCOMP = false;
446 isFULL = false;
447 isRAW = false;
448 if ( de->stwerr[se] & (1 << 16) ) isCOMP = true;
449 if ( de->stwerr[se] & (1 << 17) ) isFULL = true;
450 if ( de->stwerr[se] & (1 << 3) ) isRAW = true;
451 if ( !chdone[se] ){
452 //
453 // check for any error in the event
454 //
455 clevel2->crc[se] = 0;
456 if ( de->perror[se] == 132 ){
457 clevel2->crc[se] = 1;
458 pe++;
459 };
460 clevel2->perr[se] = 0;
461 if ( de->perror[se] != 0 ){
462 clevel2->perr[se] = 1;
463 pe++;
464 };
465 clevel2->swerr[se] = 0;
466 for (Int_t j = 0; j < 7 ; j++){
467 if ( (j != 3) && (de->stwerr[se] & (1 << j)) ){
468 clevel2->swerr[se] = 1;
469 pe++;
470 };
471 };
472 chdone[se] = 1;
473 };
474 if ( clevel2->crc[se] == 0 && (clevel1->good2 == 1 || clevel2->trigty >= 2) ){
475 pre = -1;
476 //
477 for (Int_t nn = 0; nn < 96; nn++){
478 ene[nn] = 0.;
479 dexy[l][m][nn] = de->dexy[l][m][nn] ;
480 dexyc[l][m][nn] = de->dexyc[l][m][nn] ;
481 };
482 //
483 // run over preamplifiers
484 //
485 pre = -1;
486 cbase0 = 0.;
487 for (Int_t i = 0; i < 3; i++){
488 for (Int_t j = 0; j < 2; j++){
489 pre = j + i*2;
490 //
491 // baseline check and calculation
492 //
493 if ( !isRAW ) {
494 base[l][m][pre] = de->base[l][m][pre] ;
495 cbase0 += base[l][m][pre];
496 } else {
497 //
498 // if it is a raw event and we haven't checked
499 // yet, calculate the baseline.
500 //
501 FindBaseRaw(l,m,pre);
502 cbase0 += base[l][m][pre];
503 };
504 };
505 };
506 //
507 // run over strips
508 //
509 pre = -1;
510 ener0 = 0.;
511 for (Int_t i = 0 ; i < 3 ; i++){
512 ip[i] = 0;
513 for (Int_t n = i*32 ; n < (i+1)*32 ; n++){
514 if (n%16 == 0) {
515 ck = 0;
516 done = 0;
517 doneb = 0;
518 donec = 0;
519 pre++;
520 qpre[pre] = 0.;
521 };
522 //
523 // baseline check and calculation
524 //
525 // no suitable new baseline, use old ones!
526 //
527 if ( !done ){
528 if ( (base[l][m][pre] == 31000. || base[l][m][pre] == 0.) ){
529 ck = 1;
530 if (pre%2 == 0) {
531 ip[i] = pre + 1;
532 } else {
533 ip[i] = pre - 1;
534 };
535 if ( (base[l][m][ip[i]] == 31000. || base[l][m][ip[i]] == 0.) ){
536 //
537 ck = 2;
538 if ( sbase[l][m][pre] == 31000. || sbase[l][m][pre] == 0. ) {
539 ck = 3;
540 };
541 };
542 done = 1;
543 };
544 };
545 //
546 // CALIBRATION ALGORITHM
547 //
548 if ( !doneb ){
549 switch (ck) {
550 case 0:
551 base0 = base[l][m][pre];
552 base2 = calbase[l][m][pre];
553 break;
554 case 1:
555 base0 = base[l][m][ip[i]];
556 base2 = calbase[l][m][ip[i]];
557 break;
558 case 2:
559 base0 = sbase[l][m][pre];
560 base2 = calbase[l][m][pre];
561 break;
562 case 3:
563 base0 = calbase[l][m][pre];
564 base2 = calbase[l][m][pre];
565 break;
566 };
567 base1 = calbase[l][m][pre];
568 doneb = 1;
569 };
570 ener = dexyc[l][m][n];
571 ener0 += ener;
572 clevel1->estrip[n][m][l] = 0.;
573 if ( base0>0 && base0 < 30000. ){
574 if ( !donec && (base0 - base1 + base2) != 0. ){
575 sbase[l][m][pre] = base0 - base1 + base2;
576 donec = 1;
577 };
578 if ( ener > 0. ){
579 clevel1->estrip[n][m][l] = (ener - calped[l][m][n] - base0 - base1 + base2)/mip[l][m][n] ;
580 //
581 // OK, now in estrip we have the energy deposit in MIP of all the strips for this event (at the end of loops of course)
582 //
583 qpre[pre] += clevel1->estrip[n][m][l];
584 };
585 };
586 };
587 if ( crosst ){
588 if (ck == 1){
589 if (ip[i]%2 == 0) {
590 ipre = ip[i] + 1;
591 } else {
592 ipre = ip[i] - 1;
593 };
594 for (Int_t j = ipre*16 ; j < (ipre+1)*16 ; j++){
595 clevel1->estrip[j][m][l] += (qpre[ipre] - qpre[ip[i]]) * 0.00478;
596 };
597 };
598 if (ck == 2){
599 for (Int_t j = i*32 ; j < (i+1)*32 ; j++){
600 ipre = j/16 + 1;
601 clevel1->estrip[j][m][l] += qpre[ipre] * 0.00478;
602 };
603 };
604 };
605 };
606 //
607 if ( ener0 == 0. && cbase0 == 0. && !pproblem && clevel2->perr[se] == 0){
608 if ( verbose ) printf(" L0 entry %i : calorimeter power problems! event marked as bad perr %f swerr %X view %i plane %i \n",ei,de->perror[se],de->stwerr[se],l,m);
609 pproblem = true;
610 pe++;
611 };
612 //
613 Int_t j4 = -4;
614 Int_t jjj = -3;
615 Int_t jj = -2;
616 for (Int_t j = 0 ; j < 100 ; j++){
617 jj++;
618 jjj++;
619 j4++;
620 if ( j < 96 ) ene[j] = clevel1->estrip[j][m][l];
621 if ( crosst ){
622 if ( jj >= 0 && jj < 96 ){
623 if ( jj != 0 && jj != 32 && jj != 64 ) ene[jj-1] += -clevel1->estrip[jj][m][l] * 0.01581;
624 if ( jj != 31 && jj != 63 && jj != 95 ) ene[jj+1] += -clevel1->estrip[jj][m][l] * 0.01581;
625 };
626 if ( jjj >= 0 && jjj < 96 ){
627 if ( jjj != 0 && jjj != 32 && jjj != 64 ) clevel1->estrip[jjj-1][m][l] += -ene[jjj] * 0.01581;
628 if ( jjj != 31 && jjj != 63 && jjj != 95 ) clevel1->estrip[jjj+1][m][l] += -ene[jjj] * 0.01581;
629 };
630 };
631 if ( j4 >= 0 && j4 < 96 ){
632 //
633 // NOTICE: THE FOLLOWING LINE EXCLUDE ALL STRIPS FOR WHICH THE RMS*4 IS GREATER THAN 26 !!! <=============== IMPORTANT! =================>
634 //
635 if ( obadmask[l][m][j4] == 1 || clevel1->estrip[j4][m][l] <= clevel1->emin || calrms[l][m][j4] > 26 ){
636 clevel1->estrip[j4][m][l] = 0.;
637 };
638 //
639 // code and save the energy for each strip in svstrip
640 //
641 if ( clevel1->estrip[j4][m][l] > clevel1->emin ){
642 //
643 tim = 100000.;
644 plo = m;
645 fbi = 0;
646 if ( clevel1->estrip[j4][m][l] > 0.99995 ){
647 tim = 10000.;
648 plo = m;
649 fbi = 1;
650 };
651 if ( clevel1->estrip[j4][m][l] > 9.9995 ){
652 tim = 1000.;
653 plo = 22 + m;
654 fbi = 1;
655 };
656 if ( clevel1->estrip[j4][m][l] > 99.995 ){
657 tim = 100.;
658 plo = 22 + m;
659 fbi = 0;
660 };
661 if ( clevel1->estrip[j4][m][l] > 999.95 ){
662 tim = 10.;
663 plo = 44 + m;
664 fbi = 0;
665 };
666 if ( clevel1->estrip[j4][m][l] > 9999.5 ){
667 tim = 1.;
668 plo = 66 + m;
669 fbi = 0;
670 };
671 //
672 cle = (Int_t)lroundf(tim*clevel1->estrip[j4][m][l]);
673 //
674 if ( l == 0 ){
675 //
676 // +-PPSSmmmm.mmmm
677 //
678 svstrip[istrip] = fbi*1000000000 + plo*10000000 + j4*100000 + cle;
679 } else {
680 svstrip[istrip] = -(fbi*1000000000 + plo*10000000 + j4*100000 + cle);
681 };
682 //
683 // if ( ei >= -770 ) printf(" j %i l %i m %i estrip %f \n",j4,l,m,clevel1->estrip[j4][m][l]);
684 // if ( ei >= -770 ) printf(" num lim %i fbi %i tim %f plo %i cle %i \n",numeric_limits<Int_t>::max(),fbi,tim,plo,cle);
685 // if ( ei >= -770 ) printf(" svstrip %i \n",svstrip[istrip]);
686 //
687 istrip++;
688 };
689 };
690 };
691 //
692 } else {
693 for (Int_t nn = 0; nn < 96; nn++){
694 clevel1->estrip[nn][m][l] = 0.;
695 };
696 };
697 };
698 };
699 if ( !pe ){
700 clevel2->good = 1;
701 } else {
702 clevel2->good = 0;
703 };
704 return(0);
705 }
706
707 void CaloLevel0::GetTrkVar(){
708 calol2tr();
709 }
710
711 void CaloLevel0::FillTrkVar(CaloLevel2 *ca, Int_t nutrk){
712 //
713 CaloTrkVar *t_ca = new CaloTrkVar();
714 //
715 t_ca->trkseqno = trkseqno;
716 t_ca->ncore = (Int_t)clevel2->ncore;
717 t_ca->qcore = clevel2->qcore;
718 t_ca->noint = (Int_t)clevel2->noint;
719 t_ca->ncyl = (Int_t)clevel2->ncyl;
720 t_ca->qcyl = clevel2->qcyl;
721 t_ca->qtrack = clevel2->qtrack;
722 t_ca->qtrackx = clevel2->qtrackx;
723 t_ca->qtracky = clevel2->qtracky;
724 t_ca->dxtrack = clevel2->dxtrack;
725 t_ca->dytrack = clevel2->dytrack;
726 t_ca->qlast = clevel2->qlast;
727 t_ca->nlast = (Int_t)clevel2->nlast;
728 t_ca->qpre = clevel2->qpre;
729 t_ca->npre = (Int_t)clevel2->npre;
730 t_ca->qpresh = clevel2->qpresh;
731 t_ca->npresh = (Int_t)clevel2->npresh;
732 t_ca->qtr = clevel2->qtr;
733 t_ca->ntr = (Int_t)clevel2->ntr;
734 t_ca->planetot = (Int_t)clevel2->planetot;
735 t_ca->qmean = clevel2->qmean;
736 t_ca->dX0l = clevel2->dX0l;
737 t_ca->qlow = clevel2->qlow;
738 t_ca->nlow = (Int_t)clevel2->nlow;
739 //
740 if ( trkseqno == -1 ){
741 // ca->impx = clevel2->impx;
742 // ca->impy = clevel2->impy;
743 ca->tanx[1] = clevel2->tanx;
744 ca->tany[1] = clevel2->tany;
745 ca->elen = clevel2->elen;
746 ca->selen = clevel2->selen;
747 // memcpy(ca->cibar,clevel2->cibar,sizeof(clevel2->cibar));
748 // memcpy(ca->cbar,clevel2->cbar,sizeof(clevel2->cbar));
749 memcpy(t_ca->tibar,clevel2->cibar,sizeof(clevel2->cibar));
750 memcpy(t_ca->tbar,clevel2->cbar,sizeof(clevel2->cbar));
751 memcpy(ca->planemax,clevel2->planemax,sizeof(clevel2->planemax));
752 ca->varcfit[2] = clevel2->varcfit[0];
753 ca->varcfit[3] = clevel2->varcfit[1];
754 ca->npcfit[2] = clevel2->npcfit[0];
755 ca->npcfit[3] = clevel2->npcfit[1];
756 // memcpy(ca->varcfit,clevel2->varcfit,sizeof(clevel2->varcfit));
757 // memcpy(ca->npcfit,clevel2->npcfit,sizeof(clevel2->npcfit));
758 } else {
759 memcpy(t_ca->tibar,clevel2->tibar,sizeof(clevel2->tibar));
760 memcpy(t_ca->tbar,clevel2->tbar,sizeof(clevel2->tbar));
761 };
762 //
763 //
764 if ( !(ca->CaloTrk) ) ca->CaloTrk = new TClonesArray("CaloTrkVar",1); //ELENA
765 TClonesArray &t = *ca->CaloTrk;
766 new(t[nutrk]) CaloTrkVar(*t_ca);
767 //
768 delete t_ca;
769 //
770 ClearTrkVar();
771 }
772
773 void CaloLevel0::GetCommonVar(){
774 calol2cm();
775 }
776
777 void CaloLevel0::FillCommonVar(CaloLevel1 *c1, CaloLevel2 *ca){
778 //
779 ca->good = clevel2->good;
780 if ( clevel2->trigty == 2. ){
781 ca->selftrigger = 1;
782 } else {
783 ca->selftrigger = 0;
784 };
785 //
786 ca->selftrigger += (Int_t)clevel2->wartrig;
787 //
788 memcpy(ca->perr,clevel2->perr,sizeof(clevel2->perr));
789 memcpy(ca->swerr,clevel2->swerr,sizeof(clevel2->swerr));
790 memcpy(ca->crc,clevel2->crc,sizeof(clevel2->crc));
791 ca->nstrip = (Int_t)clevel2->nstrip;
792 ca->qtot = clevel2->qtot;
793 // ca->impx = clevel2->impx;
794 // ca->impy = clevel2->impy;
795 ca->tanx[0] = clevel2->tanx;
796 ca->tany[0] = clevel2->tany;
797 ca->nx22 = (Int_t)clevel2->nx22;
798 ca->qx22 = clevel2->qx22;
799 ca->qmax = clevel2->qmax;
800 ca->elen = clevel2->elen;
801 ca->selen = clevel2->selen;
802 memcpy(ca->qq,clevel2->qq,sizeof(clevel2->qq));
803 memcpy(ca->planemax,clevel2->planemax,sizeof(clevel2->planemax));
804 ca->varcfit[0] = clevel2->varcfit[0];
805 ca->varcfit[1] = clevel2->varcfit[1];
806 ca->npcfit[0] = clevel2->npcfit[0];
807 ca->npcfit[1] = clevel2->npcfit[1];
808 ca->fitmode[0] = clevel2->fmode[0];
809 ca->fitmode[1] = clevel2->fmode[1];
810 // memcpy(ca->varcfit,clevel2->varcfit,sizeof(clevel2->varcfit));
811 // memcpy(ca->npcfit,clevel2->npcfit,sizeof(clevel2->npcfit));
812 memcpy(ca->cibar,clevel2->cibar,sizeof(clevel2->cibar));
813 memcpy(ca->cbar,clevel2->cbar,sizeof(clevel2->cbar));
814 //
815 if ( c1 ){
816 c1->istrip = istrip;
817 c1->estrip = TArrayI(istrip,svstrip);
818 };
819 //
820 }
821
822 void CaloLevel0::ClearStructs(){
823 ClearTrkVar();
824 ClearCommonVar();
825 }
826
827 void CaloLevel0::RunClose(){
828 l0tr->Delete();
829 ClearStructs();
830 //
831 memset(dexy, 0, 2*22*96*sizeof(Float_t));
832 memset(dexyc, 0, 2*22*96*sizeof(Float_t));
833 memset(base, 0, 2*22*6*sizeof(Float_t));
834 memset(sbase, 0, 2*22*6*sizeof(Float_t));
835 //
836 }
837
838 //
839 // Private methods
840 //
841
842 void CaloLevel0::ClearTrkVar(){
843 clevel2->ncore = 0;
844 clevel2->qcore = 0.;
845 clevel2->noint = 0.;
846 clevel2->ncyl = 0.;
847 clevel2->qcyl = 0.;
848 clevel2->qtrack = 0.;
849 clevel2->qtrackx = 0.;
850 clevel2->qtracky = 0.;
851 clevel2->dxtrack = 0.;
852 clevel2->dytrack = 0.;
853 clevel2->qlast = 0.;
854 clevel2->nlast = 0.;
855 clevel2->qpre = 0.;
856 clevel2->npre = 0.;
857 clevel2->qpresh = 0.;
858 clevel2->npresh = 0.;
859 clevel2->qlow = 0.;
860 clevel2->nlow = 0.;
861 clevel2->qtr = 0.;
862 clevel2->ntr = 0.;
863 clevel2->planetot = 0.;
864 clevel2->qmean = 0.;
865 clevel2->dX0l = 0.;
866 clevel2->elen = 0.;
867 clevel2->selen = 0.;
868 memset(clevel1->al_p, 0, 5*2*sizeof(Double_t));
869 memset(clevel2->tibar, 0, 2*22*sizeof(Int_t));
870 memset(clevel2->tbar, 0, 2*22*sizeof(Float_t));
871 }
872
873 void CaloLevel0::ClearCommonVar(){
874 istrip = 0;
875 clevel2->trigty = -1.;
876 clevel2->wartrig = 0.;
877 clevel2->good = 0;
878 clevel2->nstrip = 0.;
879 clevel2->qtot = 0.;
880 // clevel2->impx = 0.;
881 // clevel2->impy = 0.;
882 clevel2->tanx = 0.; // this is correct since it refers to the fortran structure
883 clevel2->tany = 0.; // this is correct since it refers to the fortran structure
884 clevel2->qmax = 0.;
885 clevel2->nx22 = 0.;
886 clevel2->qx22 = 0.;
887 memset(clevel2->perr, 0, 4*sizeof(Int_t));
888 memset(clevel2->swerr, 0, 4*sizeof(Int_t));
889 memset(clevel2->crc, 0, 4*sizeof(Int_t));
890 memset(clevel2->qq, 0, 4*sizeof(Int_t));
891 memset(clevel2->varcfit, 0, 4*sizeof(Float_t));
892 memset(clevel2->npcfit, 0, 4*sizeof(Int_t));
893 memset(clevel2->planemax, 0, 2*sizeof(Int_t));
894 memset(clevel2->fmode, 0, 2*sizeof(Int_t));
895 memset(clevel2->cibar, 0, 2*22*sizeof(Int_t));
896 memset(clevel2->cbar, 0, 2*22*sizeof(Float_t));
897 }
898
899 void CaloLevel0::ClearCalibVals(Int_t s){
900 //
901 for ( Int_t d=0 ; d<11 ;d++ ){
902 Int_t pre = -1;
903 for ( Int_t j=0; j<96 ;j++){
904 if ( j%16 == 0 ) pre++;
905 if ( s == 2 ){
906 calped[0][2*d+1][j] = 0.;
907 cstwerr[3] = 0.;
908 cperror[3] = 0.;
909 calgood[0][2*d+1][j] = 0.;
910 calthr[0][2*d+1][pre] = 0.;
911 calrms[0][2*d+1][j] = 0.;
912 calbase[0][2*d+1][pre] = 0.;
913 calvar[0][2*d+1][pre] = 0.;
914 };
915 if ( s == 3 ){
916 calped[0][2*d][j] = 0.;
917 cstwerr[1] = 0.;
918 cperror[1] = 0.;
919 calgood[0][2*d][j] = 0.;
920 calthr[0][2*d][pre] = 0.;
921 calrms[0][2*d][j] = 0.;
922 calbase[0][2*d][pre] = 0.;
923 calvar[0][2*d][pre] = 0.;
924 };
925 if ( s == 0 ){
926 calped[1][2*d][j] = 0.;
927 cstwerr[0] = 0.;
928 cperror[0] = 0.;
929 calgood[1][2*d][j] = 0.;
930 calthr[1][2*d][pre] = 0.;
931 calrms[1][2*d][j] = 0.;
932 calbase[1][2*d][pre] = 0.;
933 calvar[1][2*d][pre] = 0.;
934 };
935 if ( s == 1 ){
936 calped[1][2*d+1][j] = 0.;
937 cstwerr[2] = 0.;
938 cperror[2] = 0.;
939 calgood[1][2*d+1][j] = 0.;
940 calthr[1][2*d+1][pre] = 0.;
941 calrms[1][2*d+1][j] = 0.;
942 calbase[1][2*d+1][pre] = 0.;
943 calvar[1][2*d+1][pre] = 0.;
944 };
945 };
946 };
947 return;
948 }
949
950 Int_t CaloLevel0::Update(TSQLServer *dbc, UInt_t atime, Int_t s){
951 //
952 Int_t sgnl = 0;
953 //
954 GL_CALO_CALIB *glcalo = new GL_CALO_CALIB();
955 //
956 sgnl = 0;
957 //
958 idcalib[s] = 0;
959 fromtime[s] = 0;
960 totime[s] = 0;
961 calibno[s] = 0;
962 ClearCalibVals(s);
963 //
964 UInt_t uptime = 0;
965 //
966 sgnl = glcalo->Query_GL_CALO_CALIB(atime,uptime,s,dbc);
967 if ( sgnl < 0 ){
968 if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
969 return(sgnl);
970 };
971 //
972 idcalib[s] = glcalo->ID_ROOT_L0;
973 fromtime[s] = glcalo->FROM_TIME;
974 if ( glcalo->TO_TIME < atime ){ // calibration is corrupted and we are using the one that preceed the good one
975 totime[s] = uptime;
976 } else {
977 totime[s] = glcalo->TO_TIME;
978 };
979 // totime[s] = glcalo->TO_TIME;
980 calibno[s] = glcalo->EV_ROOT;
981 //
982 if ( totime[s] == 0 ){
983 if ( verbose ) printf(" CALORIMETER - WARNING: data with no associated calibration\n");
984 ClearCalibVals(s);
985 sgnl = 100;
986 };
987 //
988 // determine path and name and entry of the calibration file
989 //
990 GL_ROOT *glroot = new GL_ROOT();
991 if ( verbose ) printf("\n");
992 if ( verbose ) printf(" ** SECTION %i **\n",s);
993 //
994 sgnl = glroot->Query_GL_ROOT(idcalib[s],dbc);
995 if ( sgnl < 0 ){
996 if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
997 return(sgnl);
998 };
999 //
1000 stringstream name;
1001 name.str("");
1002 name << glroot->PATH.Data() << "/";
1003 name << glroot->NAME.Data();
1004 //
1005 fcalname[s] = (TString)name.str().c_str();
1006 if ( verbose ) printf(" - event at time %u. From time %u to time %u \n use file %s \n calibration at entry %i \n\n",atime,fromtime[s],totime[s],fcalname[s].Data(),calibno[s]);
1007 //
1008 sgnl = LoadCalib(s);
1009 //
1010 if ( sgnl != 0 ) return(sgnl);
1011 delete glcalo;
1012 delete glroot;
1013 //
1014 return(0);
1015 //
1016 }
1017
1018 Int_t CaloLevel0::LoadCalib(Int_t s){
1019 //
1020 ifstream myfile;
1021 myfile.open(fcalname[s].Data());
1022 if ( !myfile ){
1023 return(-107);
1024 };
1025 myfile.close();
1026 //
1027 TFile *File = new TFile(fcalname[s].Data());
1028 if ( !File ) return(-108);
1029 TTree *tr = (TTree*)File->Get("CalibCalPed");
1030 if ( !tr ) return(-109);
1031 //
1032 TBranch *calo = tr->GetBranch("CalibCalPed");
1033 //
1034 pamela::CalibCalPedEvent *ce = 0;
1035 tr->SetBranchAddress("CalibCalPed", &ce);
1036 //
1037 Long64_t ncalibs = calo->GetEntries();
1038 //
1039 if ( !ncalibs ) return(-110);
1040 //
1041 calo->GetEntry(calibno[s]);
1042 //
1043 if (ce->cstwerr[s] != 0 && ce->cperror[s] == 0 ) {
1044 for ( Int_t d=0 ; d<11 ;d++ ){
1045 Int_t pre = -1;
1046 for ( Int_t j=0; j<96 ;j++){
1047 if ( j%16 == 0 ) pre++;
1048 if ( s == 2 ){
1049 calped[0][2*d+1][j] = ce->calped[3][d][j];
1050 cstwerr[3] = ce->cstwerr[3];
1051 cperror[3] = ce->cperror[3];
1052 calgood[0][2*d+1][j] = ce->calgood[3][d][j];
1053 calthr[0][2*d+1][pre] = ce->calthr[3][d][pre];
1054 calrms[0][2*d+1][j] = ce->calrms[3][d][j];
1055 calbase[0][2*d+1][pre] = ce->calbase[3][d][pre];
1056 calvar[0][2*d+1][pre] = ce->calvar[3][d][pre];
1057 };
1058 if ( s == 3 ){
1059 calped[0][2*d][j] = ce->calped[1][d][j];
1060 cstwerr[1] = ce->cstwerr[1];
1061 cperror[1] = ce->cperror[1];
1062 calgood[0][2*d][j] = ce->calgood[1][d][j];
1063 calthr[0][2*d][pre] = ce->calthr[1][d][pre];
1064 calrms[0][2*d][j] = ce->calrms[1][d][j];
1065 calbase[0][2*d][pre] = ce->calbase[1][d][pre];
1066 calvar[0][2*d][pre] = ce->calvar[1][d][pre];
1067 };
1068 if ( s == 0 ){
1069 calped[1][2*d][j] = ce->calped[0][d][j];
1070 cstwerr[0] = ce->cstwerr[0];
1071 cperror[0] = ce->cperror[0];
1072 calgood[1][2*d][j] = ce->calgood[0][d][j];
1073 calthr[1][2*d][pre] = ce->calthr[0][d][pre];
1074 calrms[1][2*d][j] = ce->calrms[0][d][j];
1075 calbase[1][2*d][pre] = ce->calbase[0][d][pre];
1076 calvar[1][2*d][pre] = ce->calvar[0][d][pre];
1077 };
1078 if ( s == 1 ){
1079 calped[1][2*d+1][j] = ce->calped[2][d][j];
1080 cstwerr[2] = ce->cstwerr[2];
1081 cperror[2] = ce->cperror[2];
1082 calgood[1][2*d+1][j] = ce->calgood[2][d][j];
1083 calthr[1][2*d+1][pre] = ce->calthr[2][d][pre];
1084 calrms[1][2*d+1][j] = ce->calrms[2][d][j];
1085 calbase[1][2*d+1][pre] = ce->calbase[2][d][pre];
1086 calvar[1][2*d+1][pre] = ce->calvar[2][d][pre];
1087 };
1088 };
1089 };
1090 } else {
1091 if ( verbose ) printf(" CALORIMETER - ERROR: problems finding a good calibration in this file! \n\n ");
1092 return(-111);
1093 };
1094 File->Close();
1095 return(0);
1096 }

  ViewVC Help
Powered by ViewVC 1.1.23