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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.2 - (hide annotations) (download)
Fri Jun 30 09:21:54 2006 UTC (18 years, 5 months ago) by mocchiut
Branch: MAIN
CVS Tags: v1r01beta, v1r01
Changes since 1.1: +10 -0 lines
Memory leak bug fixed

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

  ViewVC Help
Powered by ViewVC 1.1.23