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

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

Parent Directory Parent Directory | Revision Log Revision Log

Revision 1.1 - (hide 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 mocchiut 1.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