/[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.10 - (hide annotations) (download)
Fri Nov 9 10:38:25 2007 UTC (17 years ago) by mocchiut
Branch: MAIN
Changes since 1.9: +10 -7 lines
Strip saturation coded in the level1 class, new method added to retrieve the information, backward compatibility preserved

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 mocchiut 1.9 memset(obadmask, 0, 2*22*96*sizeof(Int_t));
68     memset(obadpulsemask, 0, 2*22*6*sizeof(Int_t));
69     memset(ctprecor, 0, 2*22*6*sizeof(Float_t));
70     memset(ctneigcor, 0, 2*22*6*sizeof(Float_t));
71 mocchiut 1.1 calopar1 = true;
72     calopar2 = true;
73     calopar3 = true;
74     crosst = true;
75     ftcalopar1 = 0;
76     ttcalopar1 = 0;
77     ftcalopar2 = 0;
78     ttcalopar2 = 0;
79     ftcalopar3 = 0;
80     ttcalopar3 = 0;
81     }
82    
83     void CaloLevel0::SetCrossTalk(Bool_t ct){
84     crosst = ct;
85     };
86    
87 mocchiut 1.9 void CaloLevel0::SetCrossTalkType(Bool_t ct){
88     ctground = ct;
89     };
90    
91 mocchiut 1.1 void CaloLevel0::SetVerbose(Bool_t ct){
92     verbose = ct;
93     };
94    
95     /**
96     * Initialize CaloLevel0 object
97     **/
98 mocchiut 1.7 void CaloLevel0::ProcessingInit(GL_TABLES *glt, UInt_t hs, Int_t &sgnl, TTree *l0tree, Bool_t isdeb, Bool_t isverb){
99     //
100     const TString host = glt->CGetHost();
101     const TString user = glt->CGetUser();
102     const TString psw = glt->CGetPsw();
103     TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
104     if ( !dbc->IsConnected() ) throw -116;
105 mocchiut 1.8 stringstream myquery;
106     myquery.str("");
107     myquery << "SET time_zone='+0:00'";
108     dbc->Query(myquery.str().c_str());
109 mocchiut 1.1 //
110     debug = isdeb;
111     verbose = isverb;
112     //
113     l0tr=(TTree*)l0tree;
114     de = new pamela::calorimeter::CalorimeterEvent();
115     l0calo = (TBranch*)l0tr->GetBranch("Calorimeter");
116     l0tr->SetBranchAddress("Calorimeter", &de);
117     //
118     trkseqno = 0;
119     ClearStructs();
120     //
121     GL_CALO_CALIB *glcalo = new GL_CALO_CALIB();
122     //
123     sgnl = 0;
124     UInt_t uptime = 0;
125     //
126     for (Int_t s = 0; s < 4; s++){
127     idcalib[s] = 0;
128     fromtime[s] = 0;
129     totime[s] = 0;
130     calibno[s] = 0;
131     ClearCalibVals(s);
132     //
133     sgnl = glcalo->Query_GL_CALO_CALIB(hs,uptime,s,dbc);
134     if ( sgnl < 0 ){
135     if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
136     return;
137     };
138     //
139     idcalib[s] = glcalo->ID_ROOT_L0;
140     fromtime[s] = glcalo->FROM_TIME;
141     if ( glcalo->TO_TIME < hs ){ // calibration is corrupted and we are using the one that preceed the good one
142     totime[s] = uptime;
143     } else {
144     totime[s] = glcalo->TO_TIME;
145     };
146     calibno[s] = glcalo->EV_ROOT;
147     //
148     if ( totime[s] == 0 ){
149     if ( verbose ) printf(" CALORIMETER - WARNING: data with no associated calibration\n");
150     ClearCalibVals(s);
151     sgnl = 100;
152     };
153     };
154     //
155     // determine path and name and entry of the calibration file
156     //
157     GL_ROOT *glroot = new GL_ROOT();
158     if ( verbose ) printf("\n");
159     for (Int_t s = 0; s < 4; s++){
160     if ( verbose ) printf(" ** SECTION %i **\n",s);
161     if ( totime[s] > 0 ){
162     //
163     sgnl = glroot->Query_GL_ROOT(idcalib[s],dbc);
164     if ( sgnl < 0 ){
165     if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
166     return;
167     };
168     //
169     stringstream name;
170     name.str("");
171     name << glroot->PATH.Data() << "/";
172     name << glroot->NAME.Data();
173     //
174     fcalname[s] = (TString)name.str().c_str();
175     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]);
176     } else {
177     if ( verbose ) printf(" - runheader at time %u. NO CALIBRATION INCLUDE THE RUNHEADER! ",hs);
178     };
179     sgnl = LoadCalib(s);
180     if ( sgnl ) break;
181     };
182     //
183     delete glcalo;
184     delete glroot;
185 mocchiut 1.7 dbc->Close();
186     delete dbc;
187 mocchiut 1.1 //
188     return;
189     //
190     }
191    
192 mocchiut 1.7 Int_t CaloLevel0::ChkCalib(GL_TABLES *glt, UInt_t atime){
193 mocchiut 1.1 Int_t sgnl = 0;
194     for ( Int_t s = 0; s < 4; s++){
195     if ( atime > totime[s] ){
196 mocchiut 1.7 sgnl = Update(glt,atime,s);
197 mocchiut 1.1 if ( sgnl < 0 ) return(sgnl);
198     };
199     };
200     return(sgnl);
201     }
202    
203 mocchiut 1.7 Int_t CaloLevel0::ChkParam(GL_TABLES *glt, UInt_t runheader, Bool_t mechal){
204     const TString host = glt->CGetHost();
205     const TString user = glt->CGetUser();
206     const TString psw = glt->CGetPsw();
207     TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
208     if ( !dbc->IsConnected() ) throw -116;
209 mocchiut 1.8 stringstream myquery;
210     myquery.str("");
211     myquery << "SET time_zone='+0:00'";
212     dbc->Query(myquery.str().c_str());
213 mocchiut 1.7 //
214 mocchiut 1.1 stringstream calfile;
215     stringstream bmfile;
216     stringstream aligfile;
217     Int_t error = 0;
218     FILE *f = 0;
219     ifstream badfile;
220     GL_PARAM *glparam = new GL_PARAM();
221     //
222     if ( calopar1 || ( ttcalopar1 != 0 && ttcalopar1 < runheader ) ){
223 mocchiut 1.2 //
224     //
225     //
226     if ( debug ) printf(" calopar1 %i ftcalopar1 %u ttcalopar1 %u runheader %u \n",calopar1,ftcalopar1,ttcalopar1,runheader);
227     //
228 mocchiut 1.1 calopar1 = false;
229     //
230     // determine where I can find calorimeter ADC to MIP conversion file
231     //
232     if ( verbose ) printf(" Querying DB for calorimeter parameters files...\n");
233     //
234     error = 0;
235     error = glparam->Query_GL_PARAM(runheader,101,dbc);
236     if ( error < 0 ) return(error);
237     //
238     calfile.str("");
239     calfile << glparam->PATH.Data() << "/";
240     calfile << glparam->NAME.Data();
241     ftcalopar1 = glparam->FROM_TIME;
242     ttcalopar1 = glparam->TO_TIME;
243     //
244     if ( verbose ) printf("\n Using ADC to MIP conversion file: \n %s \n",calfile.str().c_str());
245     f = fopen(calfile.str().c_str(),"rb");
246     if ( !f ){
247     if ( verbose ) printf(" CALORIMETER - ERROR: no ADC to MIP file!\n");
248     return(-105);
249     };
250     //
251     for (Int_t m = 0; m < 2 ; m++ ){
252     for (Int_t k = 0; k < 22; k++ ){
253     for (Int_t l = 0; l < 96; l++ ){
254     fread(&mip[m][k][l],sizeof(mip[m][k][l]),1,f);
255 mocchiut 1.5 if ( debug ) printf(" %f \n",mip[m][k][l]);
256 mocchiut 1.1 };
257     };
258     };
259     fclose(f);
260     };
261     //
262     if ( calopar2 || ( ttcalopar2 != 0 && ttcalopar2 < runheader ) ){
263 mocchiut 1.2 if ( debug ) printf(" calopar2 %i ftcalopar2 %u ttcalopar2 %u runheader %u \n",calopar2,ftcalopar2,ttcalopar2,runheader);
264 mocchiut 1.1 calopar2 = false;
265     //
266     // determine where I can find calorimeter alignment file
267     //
268     //
269     error = 0;
270     error = glparam->Query_GL_PARAM(runheader,102,dbc);
271     if ( error < 0 ) return(error);
272     //
273     aligfile.str("");
274     aligfile << glparam->PATH.Data() << "/";
275     aligfile << glparam->NAME.Data();
276     ftcalopar2 = glparam->FROM_TIME;
277     ttcalopar2 = glparam->TO_TIME;
278     //
279 mocchiut 1.2 if ( verbose ) printf("\n Using parameter file: \n %s \n",aligfile.str().c_str());
280 mocchiut 1.1 f = fopen(aligfile.str().c_str(),"rb");
281     if ( !f ){
282 mocchiut 1.2 if ( verbose ) printf(" CALORIMETER - ERROR: no parameter file!\n");
283 mocchiut 1.1 return(-106);
284     };
285     //
286 mocchiut 1.2 if ( !mechal ){
287     //
288     fread(&clevel1->xalig,sizeof(clevel1->xalig),1,f);
289     if ( debug ) printf(" xalig = %f \n",clevel1->xalig);
290     fread(&clevel1->yalig,sizeof(clevel1->yalig),1,f);
291     if ( debug ) printf(" yalig = %f \n",clevel1->yalig);
292     fread(&clevel1->zalig,sizeof(clevel1->zalig),1,f);
293     if ( debug ) printf(" zalig = %f \n",clevel1->zalig);
294     } else {
295     if ( verbose ) printf("\n Using MECHANICAL alignement parameters \n");
296     //
297     CaloStrip cs = CaloStrip();
298     cs.UseMechanicalAlig();
299     clevel1->xalig = cs.GetXalig();
300     if ( debug ) printf(" xalig = %f \n",clevel1->xalig);
301     clevel1->yalig = cs.GetYalig();
302     if ( debug ) printf(" yalig = %f \n",clevel1->yalig);
303     clevel1->zalig = cs.GetZalig();
304     if ( debug ) printf(" zalig = %f \n",clevel1->zalig);
305     //
306     Float_t tmp = 0;
307     fread(&tmp,sizeof(clevel1->xalig),1,f);
308     fread(&tmp,sizeof(clevel1->yalig),1,f);
309     fread(&tmp,sizeof(clevel1->zalig),1,f);
310     //
311     };
312 mocchiut 1.1 fread(&clevel1->emin,sizeof(clevel1->emin),1,f);
313     if ( debug ) printf(" signal threshold = %f \n",clevel1->emin);
314     //
315     fclose(f);
316     };
317     //
318     // Load offline bad strip mask
319     //
320     if ( calopar3 || ( ttcalopar3 != 0 && ttcalopar3 < runheader ) ){
321 mocchiut 1.2 if ( debug ) printf(" calopar3 %i ftcalopar3 %u ttcalopar3 %u runheader %u \n",calopar3,ftcalopar3,ttcalopar3,runheader);
322 mocchiut 1.1 calopar3 = false;
323     //
324     // determine where I can find calorimeter alignment file
325     //
326     //
327     error = 0;
328     error = glparam->Query_GL_PARAM(runheader,103,dbc);
329     if ( error < 0 ) return(error);
330     //
331     bmfile.str("");
332     bmfile << glparam->PATH.Data() << "/";
333     bmfile << glparam->NAME.Data();
334     ftcalopar3 = glparam->FROM_TIME;
335     ttcalopar3 = glparam->TO_TIME;
336     //
337     if ( verbose ) printf("\n Using bad strip offline mask file: \n %s \n\n",bmfile.str().c_str());
338     badfile.open(bmfile.str().c_str());
339     if ( !badfile ){
340     if ( verbose ) printf(" CALORIMETER - ERROR: no bad strip offline mask file!\n");
341     return(-115);
342     };
343     //
344     Bool_t isdone = false;
345     Int_t bad = 0;
346     Int_t view = 1;
347     Int_t strip = 0;
348     Int_t plane = 21;
349     while ( !isdone ) {
350     badfile >> bad;
351     obadmask[view][plane][strip] = bad;
352     if ( debug && bad ) printf(" SETTING view %i plane %i strip %i BAD = %i \n",view,plane,strip,bad);
353     strip++;
354     if ( strip > 95 ){
355     strip = 0;
356     plane--;
357     if ( plane < 0 ){
358     plane = 21;
359     view--;
360     };
361     if ( view < 0 ) isdone = true;
362     };
363     };
364     //
365     badfile.close();
366     };
367     //
368     delete glparam;
369 mocchiut 1.7 dbc->Close();
370     delete dbc;
371 mocchiut 1.1 //
372     return(0);
373     }
374    
375 mocchiut 1.9 Int_t CaloLevel0::CalcCrossTalkCorr(GL_TABLES *glt, UInt_t runheader){
376     //
377     if ( ctground ) return(0);
378     //
379     const TString host = glt->CGetHost();
380     const TString user = glt->CGetUser();
381     const TString psw = glt->CGetPsw();
382     TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
383     if ( !dbc->IsConnected() ) throw -116;
384     stringstream myquery;
385     myquery.str("");
386     myquery << "SET time_zone='+0:00'";
387     dbc->Query(myquery.str().c_str());
388     //
389     stringstream bmfile;
390     Int_t error = 0;
391     ifstream badfile;
392     GL_PARAM *glparam = new GL_PARAM();
393     //
394     // determine where I can find file with offline bad pulser mask
395     //
396     error = 0;
397     error = glparam->Query_GL_PARAM(runheader,105,dbc);
398     if ( error < 0 ) return(error);
399     //
400     bmfile.str("");
401     bmfile << glparam->PATH.Data() << "/";
402     bmfile << glparam->NAME.Data();
403     //
404     if ( verbose ) printf("\n Using bad pulser offline mask file: \n %s \n\n",bmfile.str().c_str());
405     badfile.open(bmfile.str().c_str());
406     if ( !badfile ){
407     if ( verbose ) printf(" CALORIMETER - ERROR: no bad pulser offline mask file!\n");
408     return(-115);
409     };
410     //
411     Bool_t isdone = false;
412     Int_t bad = 0;
413     Int_t view = 1;
414     Int_t pre = 0;
415     Int_t plane = 21;
416     while ( !isdone ) {
417     badfile >> bad;
418     obadpulsemask[view][plane][pre] = bad;
419     if ( debug && bad ) printf(" SETTING view %i plane %i pre %i BAD = %i \n",view,plane,pre,bad);
420     pre++;
421     if ( pre > 5 ){
422     pre = 0;
423     plane--;
424     if ( plane < 0 ){
425     plane = 21;
426     view--;
427     };
428     if ( view < 0 ) isdone = true;
429     };
430     };
431     //
432     delete glparam;
433     badfile.close();
434     //
435     // Let's start with cross-talk correction calculation
436     //
437     GL_CALOPULSE_CALIB *glp = new GL_CALOPULSE_CALIB();
438     Float_t adcp[2][22][96];
439     Float_t adcpcal[2][22][96];
440     memset(adcp , 0, 2*22*96*sizeof(Float_t));
441     memset(adcpcal , 0, 2*22*96*sizeof(Float_t));
442     //
443     UInt_t pampli = 0;
444     for (Int_t s=0; s<4; s++){
445     //
446     // Save into matrix adcp the values of the highest pulse calibration (pulse amplitude = 2)
447     //
448     pampli = 2;
449     error = 0;
450     error = glp->Query_GL_CALOPULSE_CALIB(runheader,s,pampli,dbc);
451     if ( error < 0 ){
452     if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
453     return(error);
454     };
455     //
456     UInt_t idcalib = glp->ID_ROOT_L0;
457     UInt_t fromtime = glp->FROM_TIME;
458     UInt_t calibno = glp->EV_ROOT;
459     //
460     // determine path and name and entry of the calibration file
461     //
462     GL_ROOT *glroot = new GL_ROOT();
463     if ( verbose ) printf("\n");
464     if ( verbose ) printf(" ** SECTION %i **\n",s);
465     //
466     error = 0;
467     error = glroot->Query_GL_ROOT(idcalib,dbc);
468     if ( error < 0 ){
469     if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
470     return(error);
471     };
472     //
473     stringstream name;
474     name.str("");
475     name << glroot->PATH.Data() << "/";
476     name << glroot->NAME.Data();
477     //
478     TString fcalname = (TString)name.str().c_str();
479     ifstream myfile;
480     myfile.open(fcalname.Data());
481     if ( !myfile ){
482     return(-107);
483     };
484     myfile.close();
485     //
486     TFile *File = new TFile(fcalname.Data());
487     if ( !File ) return(-108);
488     TTree *tr = (TTree*)File->Get("CalibCalPulse2");
489     if ( !tr ) return(-119);
490     //
491     TBranch *calo = tr->GetBranch("CalibCalPulse2");
492     //
493     pamela::CalibCalPulse2Event *ce = 0;
494     tr->SetBranchAddress("CalibCalPulse2", &ce);
495     //
496     Long64_t ncalibs = calo->GetEntries();
497     //
498     if ( !ncalibs ) return(-110);
499     //
500     calo->GetEntry(calibno);
501     if ( verbose ) printf(" PULSE2 using entry %u from file %s",calibno,fcalname.Data());
502     //
503     // retrieve calibration table
504     //
505     if ( ce->pstwerr[s] && ce->pperror[s] == 0 && ce->unpackError == 0 ){
506     for ( Int_t d=0 ; d<11 ;d++ ){
507     for ( Int_t j=0; j<96 ;j++){
508     if ( s == 2 ){
509     adcp[0][2*d+1][j] = ce->calpuls[3][d][j];
510     };
511     if ( s == 3 ){
512     adcp[0][2*d][j] = ce->calpuls[1][d][j];
513     };
514     if ( s == 0 ){
515     adcp[1][2*d][j] = ce->calpuls[0][d][j];
516     };
517     if ( s == 1 ){
518     adcp[1][2*d+1][j] = ce->calpuls[2][d][j];
519     };
520     };
521     };
522     } else {
523     if ( verbose ) printf(" CALORIMETER - ERROR: problems finding a good calibration in this file! \n\n ");
524     return(-111);
525     };
526     //
527     File->Close();
528     delete glroot;
529     //
530     // Save into matrix adcpcal the calibrated values of the pulse calibration (subtraction of pulse amplitude = 0 relative to the pulse2 calibration used)
531     //
532     pampli = 0;
533     error = 0;
534     error = glp->Query_GL_CALOPULSE_CALIB(fromtime,s,pampli,dbc);
535     if ( error < 0 ){
536     if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
537     return(error);
538     };
539     //
540     idcalib = glp->ID_ROOT_L0;
541     calibno = glp->EV_ROOT;
542     //
543     // determine path and name and entry of the calibration file
544     //
545     glroot = new GL_ROOT();
546     if ( verbose ) printf("\n");
547     if ( verbose ) printf(" ** SECTION %i **\n",s);
548     //
549     error = 0;
550     error = glroot->Query_GL_ROOT(idcalib,dbc);
551     if ( error < 0 ){
552     if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
553     return(error);
554     };
555     //
556     name.str("");
557     name << glroot->PATH.Data() << "/";
558     name << glroot->NAME.Data();
559     //
560     fcalname = (TString)name.str().c_str();
561     myfile.open(fcalname.Data());
562     if ( !myfile ){
563     return(-107);
564     };
565     myfile.close();
566     //
567     TFile *File1 = new TFile(fcalname.Data());
568     if ( !File1 ) return(-108);
569     TTree *tr1 = (TTree*)File1->Get("CalibCalPulse1");
570     if ( !tr1 ) return(-120);
571     //
572     TBranch *calo1 = tr1->GetBranch("CalibCalPulse1");
573     //
574     pamela::CalibCalPulse1Event *ce1 = 0;
575     tr1->SetBranchAddress("CalibCalPulse1", &ce1);
576     //
577     ncalibs = calo1->GetEntries();
578     //
579     if ( !ncalibs ) return(-110);
580     //
581     calo1->GetEntry(calibno);
582     if ( verbose ) printf(" PULSE1 using entry %u from file %s",calibno,fcalname.Data());
583     //
584     // retrieve calibration table
585     //
586     if ( ce1->pstwerr[s] && ce1->pperror[s] == 0 && ce1->unpackError == 0 ){
587     for ( Int_t d=0 ; d<11 ;d++ ){
588     for ( Int_t j=0; j<96 ;j++){
589     if ( s == 2 ){
590     adcpcal[0][2*d+1][j] = adcp[0][2*d+1][j] - ce1->calpuls[3][d][j];
591     };
592     if ( s == 3 ){
593     adcpcal[0][2*d][j] = adcp[0][2*d][j] - ce1->calpuls[1][d][j];
594     };
595     if ( s == 0 ){
596     adcpcal[1][2*d][j] = adcp[1][2*d][j] - ce1->calpuls[0][d][j];
597     };
598     if ( s == 1 ){
599     adcpcal[1][2*d+1][j] = adcp[1][2*d+1][j] - ce1->calpuls[2][d][j];
600     };
601     };
602     };
603     } else {
604     if ( verbose ) printf(" CALORIMETER - ERROR: problems finding a good calibration in this file! \n\n ");
605     return(-111);
606     };
607     //
608     File1->Close();
609     //
610     delete glroot;
611     //
612     };// loop on the four sections
613     //
614     //
615     delete glp;
616     dbc->Close();
617     delete dbc;
618     //
619     // Ok, now we can try to calculate the cross-talk correction for each pre-amplifier
620     //
621     for ( Int_t v=0; v<2; v++){
622     if ( debug ) printf(" \n\n NEW VIEW \n");
623     for ( Int_t p=0; p<22; p++){
624     for ( Int_t npre=0; npre<6; npre++){
625     ctprecor[v][p][npre] = 1000.;
626     ctneigcor[v][p][npre] = 1000.;
627     Int_t str0=npre*16;
628     Int_t str16= -1 + (1+npre)*16;
629     //
630     UInt_t neigc = 0;
631     UInt_t farc = 0;
632     UInt_t pulsc = 0;
633     Float_t sigpulsed = 0.;
634     Float_t neigbase = 0.;
635     Float_t farbase = 0.;
636     //
637     // Loop over the strip of the pre and sum all signal far away from pulsed strip, signal in the neighbour(s) strip(s) and save the pulsed signal
638     // moreover count the number of strips in each case
639     //
640     for (Int_t s=str0; s<=str16; s++){
641     if ( adcpcal[v][p][s] > 10000.){
642     sigpulsed = adcpcal[v][p][s];
643     pulsc++;
644     if ( s > str0 ){
645     neigbase += adcpcal[v][p][s-1];
646     neigc++;
647     farbase -= adcpcal[v][p][s-1];
648     farc--;
649     };
650     if ( s < str16 ){
651     neigbase += adcpcal[v][p][s+1];
652     neigc++;
653     farbase -= adcpcal[v][p][s+1];
654     farc--;
655     };
656     } else {
657     farc++;
658     farbase += adcpcal[v][p][s];
659     };
660     };
661     //
662     // Now calculate the corrections
663     //
664     Float_t avefarbase = 0.;
665     if ( farc ) avefarbase = farbase/(Float_t)farc;
666     Float_t aveneigbase = 0.;
667     if ( neigc ) aveneigbase = neigbase/(Float_t)neigc;
668     //
669     if ( pulsc == 1 && farc && neigc ){
670     ctprecor[v][p][npre] = -avefarbase/(sigpulsed+fabs(avefarbase));
671     ctneigcor[v][p][npre] = fabs(aveneigbase-avefarbase)/(sigpulsed+fabs(avefarbase));
672     if ( debug ) printf(" Cross-talk correction View %i Plane %i Pre %i : pre-correction: %f neighbour strips correction %f \n",v,p,npre,ctprecor[v][p][npre],ctneigcor[v][p][npre]);
673     } else {
674     //
675     // did not find the pulsed strip or more than one pulsed strip found!
676     //
677     if ( debug ) printf(" Problems finding the cross-talk corrections: \n View %i Plane %i Pre %i number of pulsed strip %i \n Average faraway baseline %f number of strips %i Average neighbour baseline %f number of neighbour strips %i \n",v,p,npre,pulsc,avefarbase,farc,aveneigbase,neigc);
678     //
679     };
680     };
681     if ( debug ) printf(" \n ==================== \n");
682     };
683     };
684     //
685     // Check the calculated corrections
686     //
687     Int_t opre=0;
688     Int_t ppre=0;
689     Bool_t found = false;
690     for ( Int_t v=0; v<2; v++){
691     for ( Int_t p=0; p<22; p++){
692     for ( Int_t npre=0; npre<6; npre++){
693     if ( ctprecor[v][p][npre] == 1000. || ctneigcor[v][p][npre] == 1000. || obadpulsemask[v][p][npre] != 0 ){
694     if ( debug ) printf(" Cross-talk correction CHANGED for view %i Plane %i Pre %i\n BEFORE: pre-correction: %f neighbour strips correction %f \n",v,p,npre,ctprecor[v][p][npre],ctneigcor[v][p][npre]);
695     if ( npre%2 ){
696     opre = npre-1;
697     } else {
698     opre = npre+1;
699     };
700     if ( ctprecor[v][p][opre] == 1000. || ctneigcor[v][p][opre] == 1000. || obadpulsemask[v][p][opre] != 0 ){
701     ppre=0;
702     found = false;
703     while ( ppre < 6 ){
704     if ( ctprecor[v][p][ppre] != 1000. && ctneigcor[v][p][ppre] != 1000. && !obadpulsemask[v][p][ppre] ){
705     found = true;
706     ctprecor[v][p][npre] = ctprecor[v][p][ppre];
707     ctneigcor[v][p][npre] = ctneigcor[v][p][ppre];
708     break;
709     };
710     ppre++;
711     };
712     if ( !found ){
713     if ( verbose ) printf(" WARNING: cannot find a good cross-talk correction for view %i plane %i pre %i \n Setting to default values 0.002 0.002\n",v,p,npre);
714     ctprecor[v][p][npre] = 0.002;
715     ctneigcor[v][p][npre] = 0.002;
716     };
717     } else {
718     ctprecor[v][p][npre] = ctprecor[v][p][opre];
719     ctneigcor[v][p][npre] = ctneigcor[v][p][opre];
720     };
721     if ( debug ) printf(" AFTER: pre-correction: %f neighbour strips correction %f \n",ctprecor[v][p][npre],ctneigcor[v][p][npre]);
722     };
723     };
724     };
725     };
726     //
727     return(0);
728     };
729 mocchiut 1.1
730     void CaloLevel0::FindBaseRaw(Int_t l, Int_t m, Int_t pre){
731     Float_t minstrip = 100000.;
732     Float_t rms = 0.;
733     base[l][m][pre] = 0.;
734     for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
735     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.) {
736     minstrip = dexy[l][m][e]-calped[l][m][e];
737     rms = calthr[l][m][pre];
738     };
739     };
740     if ( minstrip != 100000. ) {
741     Float_t strip6s = 0.;
742     for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
743     if ( (dexy[l][m][e]-calped[l][m][e]) >= minstrip && (dexy[l][m][e]-calped[l][m][e]) <= (minstrip+rms) ) {
744     strip6s += 1.;
745     base[l][m][pre] += (dexy[l][m][e] - calped[l][m][e]);
746     };
747     //
748     // compression
749     //
750     if ( abs((int)(dexy[l][m][e]-calped[l][m][e])) <= (minstrip+rms) ) {
751     dexyc[l][m][e] = 0.;
752     } else {
753     dexyc[l][m][e] = dexy[l][m][e];
754     };
755     };
756     if ( strip6s >= 9. ){
757     Double_t arro = base[l][m][pre]/strip6s;
758     Float_t deci = 1000.*((float)arro - float(int(arro)));
759     if ( deci < 500. ) {
760     arro = double(int(arro));
761     } else {
762     arro = 1. + double(int(arro));
763     };
764     base[l][m][pre] = arro;
765     } else {
766     base[l][m][pre] = 31000.;
767     for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
768     dexyc[l][m][e] = dexy[l][m][e];
769     };
770     };
771     } else {
772     base[l][m][pre] = 31000.;
773     };
774     }
775    
776     Int_t CaloLevel0::Calibrate(Int_t ei){
777     //
778     // get entry ei
779     //
780     l0calo->GetEntry(ei);
781     //
782     // if it was not a selftrigger event, could it ever been a selftrigger event? if so trigty = 3.
783     //
784     Int_t val = 0;
785     Int_t del = 1100;
786 mocchiut 1.6 for (Int_t sec = 0; sec < 4; sec++){
787     for (Int_t dsec = 0; dsec < 7; dsec++){
788     val = (Int_t)de->calselftrig[sec][dsec];
789     del = delay(val);
790     clevel2->selfdelay[sec][dsec] = del;
791     };
792     };
793     val = 0;
794     del = 1100;
795 mocchiut 1.1 if ( clevel2->trigty != 2. ){
796     Bool_t ck = false;
797     for (Int_t sec = 0; sec < 4; sec++){
798     val = (Int_t)de->calselftrig[sec][6];
799     del = delay(val);
800     if ( del < 1100 ){
801     clevel2->wartrig = 0.;
802     clevel2->trigty = 3.;
803     ck = true;
804     break;
805     };
806     };
807     if ( !ck ) clevel2->wartrig = 100.;
808     } else {
809     Bool_t ck = false;
810     for (Int_t sec = 0; sec < 4; sec++){
811     val = (Int_t)de->calselftrig[sec][6];
812     del = delay(val);
813     if ( del < 1100 ){
814     clevel2->wartrig = 0.;
815     ck = true;
816     };
817     };
818     if ( !ck ) clevel2->wartrig = 100.;
819     };
820     //
821     Int_t se = 5;
822     Int_t done = 0;
823     Int_t pre = -1;
824     Bool_t isCOMP = false;
825     Bool_t isFULL = false;
826     Bool_t isRAW = false;
827     Float_t ener;
828     Int_t doneb = 0;
829     Int_t donec = 0;
830     Int_t ck = 0;
831     Int_t ipre = 0;
832     Int_t ip[3] = {0};
833     Float_t base0, base1, base2;
834     base0 = 0.;
835     base1 = 0.;
836     base2 = 0.;
837     Float_t qpre[6] = {0.,0.,0.,0.,0.,0.};
838     Float_t ene[96];
839     Int_t chdone[4] = {0,0,0,0};
840     Int_t pe = 0;
841     //
842     Float_t ener0 = 0.;
843     Float_t cbase0 = 0.;
844     Bool_t pproblem = false;
845     //
846     Float_t tim = 0.;
847     Int_t plo = 0;
848     Int_t fbi = 0;
849     Int_t cle = 0;
850     //
851     // run over views and planes
852     //
853     for (Int_t l = 0; l < 2; l++){
854     for (Int_t m = 0; m < 22; m++){
855     //
856 mocchiut 1.2 // determine the section number
857 mocchiut 1.1 //
858     se = 5;
859 mocchiut 1.2 if (l == 0 && m%2 == 0) se = 3;
860 mocchiut 1.1 if (l == 0 && m%2 != 0) se = 2;
861 mocchiut 1.2 if (l == 1 && m%2 != 0) se = 1;
862     if (l == 1 && m%2 == 0) se = 0;
863 mocchiut 1.1 //
864     // determine what kind of event we are going to analyze
865     //
866     isCOMP = false;
867     isFULL = false;
868     isRAW = false;
869     if ( de->stwerr[se] & (1 << 16) ) isCOMP = true;
870     if ( de->stwerr[se] & (1 << 17) ) isFULL = true;
871     if ( de->stwerr[se] & (1 << 3) ) isRAW = true;
872     if ( !chdone[se] ){
873     //
874     // check for any error in the event
875     //
876     clevel2->crc[se] = 0;
877     if ( de->perror[se] == 132 ){
878     clevel2->crc[se] = 1;
879     pe++;
880     };
881     clevel2->perr[se] = 0;
882     if ( de->perror[se] != 0 ){
883     clevel2->perr[se] = 1;
884     pe++;
885     };
886     clevel2->swerr[se] = 0;
887     for (Int_t j = 0; j < 7 ; j++){
888     if ( (j != 3) && (de->stwerr[se] & (1 << j)) ){
889     clevel2->swerr[se] = 1;
890     pe++;
891     };
892     };
893     chdone[se] = 1;
894     };
895     if ( clevel2->crc[se] == 0 && (clevel1->good2 == 1 || clevel2->trigty >= 2) ){
896     pre = -1;
897     //
898     for (Int_t nn = 0; nn < 96; nn++){
899     ene[nn] = 0.;
900     dexy[l][m][nn] = de->dexy[l][m][nn] ;
901     dexyc[l][m][nn] = de->dexyc[l][m][nn] ;
902     };
903     //
904     // run over preamplifiers
905     //
906     pre = -1;
907     cbase0 = 0.;
908     for (Int_t i = 0; i < 3; i++){
909     for (Int_t j = 0; j < 2; j++){
910     pre = j + i*2;
911     //
912     // baseline check and calculation
913     //
914     if ( !isRAW ) {
915     base[l][m][pre] = de->base[l][m][pre] ;
916     cbase0 += base[l][m][pre];
917     } else {
918     //
919     // if it is a raw event and we haven't checked
920     // yet, calculate the baseline.
921     //
922     FindBaseRaw(l,m,pre);
923     cbase0 += base[l][m][pre];
924     };
925     };
926     };
927     //
928     // run over strips
929     //
930     pre = -1;
931     ener0 = 0.;
932     for (Int_t i = 0 ; i < 3 ; i++){
933     ip[i] = 0;
934     for (Int_t n = i*32 ; n < (i+1)*32 ; n++){
935     if (n%16 == 0) {
936     ck = 0;
937     done = 0;
938     doneb = 0;
939     donec = 0;
940     pre++;
941     qpre[pre] = 0.;
942     };
943     //
944     // baseline check and calculation
945     //
946     // no suitable new baseline, use old ones!
947     //
948     if ( !done ){
949     if ( (base[l][m][pre] == 31000. || base[l][m][pre] == 0.) ){
950     ck = 1;
951     if (pre%2 == 0) {
952     ip[i] = pre + 1;
953     } else {
954     ip[i] = pre - 1;
955     };
956 mocchiut 1.3 if ( (base[l][m][ip[i]] == 31000. || base[l][m][ip[i]] == 0. || !crosst ) ){
957 mocchiut 1.1 //
958     ck = 2;
959     if ( sbase[l][m][pre] == 31000. || sbase[l][m][pre] == 0. ) {
960     ck = 3;
961     };
962     };
963     done = 1;
964     };
965     };
966     //
967     // CALIBRATION ALGORITHM
968     //
969     if ( !doneb ){
970 mocchiut 1.3 if ( debug ) printf(" ck is %i \n",ck);
971 mocchiut 1.1 switch (ck) {
972     case 0:
973     base0 = base[l][m][pre];
974     base2 = calbase[l][m][pre];
975 mocchiut 1.3 if ( debug ) printf(" base0 = base l m pre = %f base2 = calbase l m pre = %f \n",base[l][m][pre],calbase[l][m][pre]);
976 mocchiut 1.1 break;
977     case 1:
978     base0 = base[l][m][ip[i]];
979     base2 = calbase[l][m][ip[i]];
980 mocchiut 1.3 if ( debug ) printf(" base0 = base l m ip(i) = %f base2 = calbase l m ip(i) = %f \n",base[l][m][ip[i]],calbase[l][m][ip[i]]);
981 mocchiut 1.1 break;
982     case 2:
983     base0 = sbase[l][m][pre];
984 mocchiut 1.3 base2 = calbase[l][m][pre];
985     if ( debug ) printf(" base0 = sbase l m pre = %f base2 = calbase l m pre = %f \n",sbase[l][m][pre],calbase[l][m][pre]);
986 mocchiut 1.1 break;
987     case 3:
988     base0 = calbase[l][m][pre];
989     base2 = calbase[l][m][pre];
990 mocchiut 1.3 if ( debug ) printf(" base0 = calbase l m pre = %f base2 = calbase l m pre = %f \n",calbase[l][m][pre],calbase[l][m][pre]);
991 mocchiut 1.1 break;
992     };
993     base1 = calbase[l][m][pre];
994     doneb = 1;
995     };
996     ener = dexyc[l][m][n];
997     ener0 += ener;
998     clevel1->estrip[n][m][l] = 0.;
999     if ( base0>0 && base0 < 30000. ){
1000 mocchiut 1.3 // if ( !donec && (base0 - base1 + base2) != 0. ){
1001     // sbase[l][m][pre] = base0 - base1 + base2;
1002     if ( !donec && (base0 + base1 - base2) != 0. ){
1003     sbase[l][m][pre] = base0 + base1 - base2;
1004 mocchiut 1.1 donec = 1;
1005     };
1006     if ( ener > 0. ){
1007     clevel1->estrip[n][m][l] = (ener - calped[l][m][n] - base0 - base1 + base2)/mip[l][m][n] ;
1008     //
1009     // 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)
1010     //
1011     qpre[pre] += clevel1->estrip[n][m][l];
1012 mocchiut 1.3 //
1013     //
1014 mocchiut 1.1 };
1015     };
1016     };
1017     if ( crosst ){
1018     if (ck == 1){
1019     if (ip[i]%2 == 0) {
1020     ipre = ip[i] + 1;
1021     } else {
1022     ipre = ip[i] - 1;
1023     };
1024     for (Int_t j = ipre*16 ; j < (ipre+1)*16 ; j++){
1025 mocchiut 1.9 if ( !ctground ){
1026     clevel1->estrip[j][m][l] += (qpre[ipre] - qpre[ip[i]]) * ctprecor[l][m][ipre];
1027     } else {
1028     clevel1->estrip[j][m][l] += (qpre[ipre] - qpre[ip[i]]) * 0.00478;
1029     };
1030 mocchiut 1.1 };
1031     };
1032     if (ck == 2){
1033     for (Int_t j = i*32 ; j < (i+1)*32 ; j++){
1034     ipre = j/16 + 1;
1035 mocchiut 1.9 if ( !ctground ){
1036     clevel1->estrip[j][m][l] += qpre[ipre] * ctprecor[l][m][ipre];
1037     } else {
1038     clevel1->estrip[j][m][l] += qpre[ipre] * 0.00478;
1039     };
1040 mocchiut 1.1 };
1041     };
1042     };
1043     };
1044     //
1045     if ( ener0 == 0. && cbase0 == 0. && !pproblem && clevel2->perr[se] == 0){
1046     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);
1047     pproblem = true;
1048     pe++;
1049     };
1050     //
1051     Int_t j4 = -4;
1052     Int_t jjj = -3;
1053     Int_t jj = -2;
1054 mocchiut 1.9 Int_t jjpre = -1;
1055     Int_t jjjpre = -1;
1056     for (Int_t j = 0 ; j < 100 ; j++){
1057 mocchiut 1.1 jj++;
1058     jjj++;
1059     j4++;
1060     if ( j < 96 ) ene[j] = clevel1->estrip[j][m][l];
1061     if ( crosst ){
1062     if ( jj >= 0 && jj < 96 ){
1063 mocchiut 1.9 if ( !ctground ){
1064     if ( jj%16 == 0 ) jjpre++;
1065     if ( jj != 0 && jj != 32 && jj != 64 ) ene[jj-1] += -clevel1->estrip[jj][m][l] * ctneigcor[l][m][jjpre];
1066     if ( jj != 31 && jj != 63 && jj != 95 ) ene[jj+1] += -clevel1->estrip[jj][m][l] * ctneigcor[l][m][jjpre];
1067     } else {
1068     if ( jj != 0 && jj != 32 && jj != 64 ) ene[jj-1] += -clevel1->estrip[jj][m][l] * 0.01581;
1069     if ( jj != 31 && jj != 63 && jj != 95 ) ene[jj+1] += -clevel1->estrip[jj][m][l] * 0.01581;
1070     };
1071 mocchiut 1.1 };
1072     if ( jjj >= 0 && jjj < 96 ){
1073 mocchiut 1.9 if ( !ctground ){
1074     if ( jjj%16 == 0 ) jjjpre++;
1075     if ( jjj != 0 && jjj != 32 && jjj != 64 ) clevel1->estrip[jjj-1][m][l] += -ene[jjj] * ctneigcor[l][m][jjjpre];
1076     if ( jjj != 31 && jjj != 63 && jjj != 95 ) clevel1->estrip[jjj+1][m][l] += -ene[jjj] * ctneigcor[l][m][jjjpre];
1077     } else {
1078     if ( jjj != 0 && jjj != 32 && jjj != 64 ) clevel1->estrip[jjj-1][m][l] += -ene[jjj] * 0.01581;
1079     if ( jjj != 31 && jjj != 63 && jjj != 95 ) clevel1->estrip[jjj+1][m][l] += -ene[jjj] * 0.01581;
1080     };
1081 mocchiut 1.1 };
1082     };
1083     if ( j4 >= 0 && j4 < 96 ){
1084     //
1085     // NOTICE: THE FOLLOWING LINE EXCLUDE ALL STRIPS FOR WHICH THE RMS*4 IS GREATER THAN 26 !!! <=============== IMPORTANT! =================>
1086     //
1087     if ( obadmask[l][m][j4] == 1 || clevel1->estrip[j4][m][l] <= clevel1->emin || calrms[l][m][j4] > 26 ){
1088     clevel1->estrip[j4][m][l] = 0.;
1089     };
1090 mocchiut 1.10 //
1091 mocchiut 1.1 // code and save the energy for each strip in svstrip
1092     //
1093     if ( clevel1->estrip[j4][m][l] > clevel1->emin ){
1094     //
1095 mocchiut 1.10 Float_t savel1 = clevel1->estrip[j4][m][l];
1096     if ( dexyc[j4][m][l] == 32767. ) savel1 += 5000.;
1097     //
1098 mocchiut 1.1 tim = 100000.;
1099     plo = m;
1100     fbi = 0;
1101 mocchiut 1.10 if ( savel1 > 0.99995 ){
1102 mocchiut 1.1 tim = 10000.;
1103     plo = m;
1104     fbi = 1;
1105     };
1106 mocchiut 1.10 if ( savel1 > 9.9995 ){
1107 mocchiut 1.1 tim = 1000.;
1108     plo = 22 + m;
1109     fbi = 1;
1110     };
1111 mocchiut 1.10 if ( savel1 > 99.995 ){
1112 mocchiut 1.1 tim = 100.;
1113     plo = 22 + m;
1114     fbi = 0;
1115     };
1116 mocchiut 1.10 if ( savel1 > 999.95 ){
1117 mocchiut 1.1 tim = 10.;
1118     plo = 44 + m;
1119     fbi = 0;
1120     };
1121 mocchiut 1.10 if ( savel1 > 9999.5 ){
1122 mocchiut 1.1 tim = 1.;
1123     plo = 66 + m;
1124     fbi = 0;
1125     };
1126     //
1127 mocchiut 1.10 cle = (Int_t)lroundf(tim*savel1);
1128 mocchiut 1.1 //
1129     if ( l == 0 ){
1130     //
1131     // +-PPSSmmmm.mmmm
1132     //
1133     svstrip[istrip] = fbi*1000000000 + plo*10000000 + j4*100000 + cle;
1134     } else {
1135     svstrip[istrip] = -(fbi*1000000000 + plo*10000000 + j4*100000 + cle);
1136     };
1137     //
1138     // if ( ei >= -770 ) printf(" j %i l %i m %i estrip %f \n",j4,l,m,clevel1->estrip[j4][m][l]);
1139     // 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);
1140     // if ( ei >= -770 ) printf(" svstrip %i \n",svstrip[istrip]);
1141     //
1142     istrip++;
1143     };
1144     };
1145     };
1146     //
1147     } else {
1148     for (Int_t nn = 0; nn < 96; nn++){
1149     clevel1->estrip[nn][m][l] = 0.;
1150     };
1151     };
1152     };
1153     };
1154     if ( !pe ){
1155     clevel2->good = 1;
1156     } else {
1157     clevel2->good = 0;
1158     };
1159     return(0);
1160     }
1161    
1162     void CaloLevel0::GetTrkVar(){
1163     calol2tr();
1164     }
1165    
1166     void CaloLevel0::FillTrkVar(CaloLevel2 *ca, Int_t nutrk){
1167     //
1168     CaloTrkVar *t_ca = new CaloTrkVar();
1169     //
1170     t_ca->trkseqno = trkseqno;
1171     t_ca->ncore = (Int_t)clevel2->ncore;
1172     t_ca->qcore = clevel2->qcore;
1173     t_ca->noint = (Int_t)clevel2->noint;
1174     t_ca->ncyl = (Int_t)clevel2->ncyl;
1175     t_ca->qcyl = clevel2->qcyl;
1176     t_ca->qtrack = clevel2->qtrack;
1177     t_ca->qtrackx = clevel2->qtrackx;
1178     t_ca->qtracky = clevel2->qtracky;
1179     t_ca->dxtrack = clevel2->dxtrack;
1180     t_ca->dytrack = clevel2->dytrack;
1181     t_ca->qlast = clevel2->qlast;
1182     t_ca->nlast = (Int_t)clevel2->nlast;
1183     t_ca->qpre = clevel2->qpre;
1184     t_ca->npre = (Int_t)clevel2->npre;
1185     t_ca->qpresh = clevel2->qpresh;
1186     t_ca->npresh = (Int_t)clevel2->npresh;
1187     t_ca->qtr = clevel2->qtr;
1188     t_ca->ntr = (Int_t)clevel2->ntr;
1189     t_ca->planetot = (Int_t)clevel2->planetot;
1190     t_ca->qmean = clevel2->qmean;
1191     t_ca->dX0l = clevel2->dX0l;
1192     t_ca->qlow = clevel2->qlow;
1193     t_ca->nlow = (Int_t)clevel2->nlow;
1194     //
1195     if ( trkseqno == -1 ){
1196     // ca->impx = clevel2->impx;
1197     // ca->impy = clevel2->impy;
1198     ca->tanx[1] = clevel2->tanx;
1199     ca->tany[1] = clevel2->tany;
1200     ca->elen = clevel2->elen;
1201     ca->selen = clevel2->selen;
1202     // memcpy(ca->cibar,clevel2->cibar,sizeof(clevel2->cibar));
1203     // memcpy(ca->cbar,clevel2->cbar,sizeof(clevel2->cbar));
1204     memcpy(t_ca->tibar,clevel2->cibar,sizeof(clevel2->cibar));
1205     memcpy(t_ca->tbar,clevel2->cbar,sizeof(clevel2->cbar));
1206     memcpy(ca->planemax,clevel2->planemax,sizeof(clevel2->planemax));
1207 mocchiut 1.6 memcpy(ca->selfdelay,clevel2->selfdelay,sizeof(clevel2->selfdelay));
1208 mocchiut 1.1 ca->varcfit[2] = clevel2->varcfit[0];
1209     ca->varcfit[3] = clevel2->varcfit[1];
1210     ca->npcfit[2] = clevel2->npcfit[0];
1211     ca->npcfit[3] = clevel2->npcfit[1];
1212     // memcpy(ca->varcfit,clevel2->varcfit,sizeof(clevel2->varcfit));
1213     // memcpy(ca->npcfit,clevel2->npcfit,sizeof(clevel2->npcfit));
1214     } else {
1215     memcpy(t_ca->tibar,clevel2->tibar,sizeof(clevel2->tibar));
1216     memcpy(t_ca->tbar,clevel2->tbar,sizeof(clevel2->tbar));
1217     };
1218     //
1219     //
1220     if ( !(ca->CaloTrk) ) ca->CaloTrk = new TClonesArray("CaloTrkVar",1); //ELENA
1221     TClonesArray &t = *ca->CaloTrk;
1222     new(t[nutrk]) CaloTrkVar(*t_ca);
1223     //
1224     delete t_ca;
1225     //
1226     ClearTrkVar();
1227     }
1228    
1229     void CaloLevel0::GetCommonVar(){
1230     calol2cm();
1231     }
1232    
1233     void CaloLevel0::FillCommonVar(CaloLevel1 *c1, CaloLevel2 *ca){
1234     //
1235     ca->good = clevel2->good;
1236     if ( clevel2->trigty == 2. ){
1237     ca->selftrigger = 1;
1238     } else {
1239     ca->selftrigger = 0;
1240     };
1241     //
1242     ca->selftrigger += (Int_t)clevel2->wartrig;
1243     //
1244     memcpy(ca->perr,clevel2->perr,sizeof(clevel2->perr));
1245     memcpy(ca->swerr,clevel2->swerr,sizeof(clevel2->swerr));
1246     memcpy(ca->crc,clevel2->crc,sizeof(clevel2->crc));
1247     ca->nstrip = (Int_t)clevel2->nstrip;
1248     ca->qtot = clevel2->qtot;
1249     // ca->impx = clevel2->impx;
1250     // ca->impy = clevel2->impy;
1251     ca->tanx[0] = clevel2->tanx;
1252     ca->tany[0] = clevel2->tany;
1253     ca->nx22 = (Int_t)clevel2->nx22;
1254     ca->qx22 = clevel2->qx22;
1255     ca->qmax = clevel2->qmax;
1256     ca->elen = clevel2->elen;
1257     ca->selen = clevel2->selen;
1258     memcpy(ca->qq,clevel2->qq,sizeof(clevel2->qq));
1259     memcpy(ca->planemax,clevel2->planemax,sizeof(clevel2->planemax));
1260 mocchiut 1.6 memcpy(ca->selfdelay,clevel2->selfdelay,sizeof(clevel2->selfdelay));
1261 mocchiut 1.1 ca->varcfit[0] = clevel2->varcfit[0];
1262     ca->varcfit[1] = clevel2->varcfit[1];
1263     ca->npcfit[0] = clevel2->npcfit[0];
1264     ca->npcfit[1] = clevel2->npcfit[1];
1265     ca->fitmode[0] = clevel2->fmode[0];
1266     ca->fitmode[1] = clevel2->fmode[1];
1267     // memcpy(ca->varcfit,clevel2->varcfit,sizeof(clevel2->varcfit));
1268     // memcpy(ca->npcfit,clevel2->npcfit,sizeof(clevel2->npcfit));
1269     memcpy(ca->cibar,clevel2->cibar,sizeof(clevel2->cibar));
1270     memcpy(ca->cbar,clevel2->cbar,sizeof(clevel2->cbar));
1271     //
1272     if ( c1 ){
1273     c1->istrip = istrip;
1274     c1->estrip = TArrayI(istrip,svstrip);
1275     };
1276     //
1277     }
1278    
1279     void CaloLevel0::ClearStructs(){
1280     ClearTrkVar();
1281     ClearCommonVar();
1282     }
1283    
1284     void CaloLevel0::RunClose(){
1285     l0tr->Delete();
1286     ClearStructs();
1287     //
1288     memset(dexy, 0, 2*22*96*sizeof(Float_t));
1289     memset(dexyc, 0, 2*22*96*sizeof(Float_t));
1290     memset(base, 0, 2*22*6*sizeof(Float_t));
1291     memset(sbase, 0, 2*22*6*sizeof(Float_t));
1292 mocchiut 1.9 memset(ctprecor, 0, 2*22*6*sizeof(Float_t));
1293     memset(ctneigcor, 0, 2*22*6*sizeof(Float_t));
1294 mocchiut 1.1 //
1295     }
1296    
1297     //
1298     // Private methods
1299     //
1300    
1301     void CaloLevel0::ClearTrkVar(){
1302     clevel2->ncore = 0;
1303     clevel2->qcore = 0.;
1304     clevel2->noint = 0.;
1305     clevel2->ncyl = 0.;
1306     clevel2->qcyl = 0.;
1307     clevel2->qtrack = 0.;
1308     clevel2->qtrackx = 0.;
1309     clevel2->qtracky = 0.;
1310     clevel2->dxtrack = 0.;
1311     clevel2->dytrack = 0.;
1312     clevel2->qlast = 0.;
1313     clevel2->nlast = 0.;
1314     clevel2->qpre = 0.;
1315     clevel2->npre = 0.;
1316     clevel2->qpresh = 0.;
1317     clevel2->npresh = 0.;
1318     clevel2->qlow = 0.;
1319     clevel2->nlow = 0.;
1320     clevel2->qtr = 0.;
1321     clevel2->ntr = 0.;
1322     clevel2->planetot = 0.;
1323     clevel2->qmean = 0.;
1324     clevel2->dX0l = 0.;
1325     clevel2->elen = 0.;
1326     clevel2->selen = 0.;
1327     memset(clevel1->al_p, 0, 5*2*sizeof(Double_t));
1328     memset(clevel2->tibar, 0, 2*22*sizeof(Int_t));
1329     memset(clevel2->tbar, 0, 2*22*sizeof(Float_t));
1330     }
1331    
1332     void CaloLevel0::ClearCommonVar(){
1333     istrip = 0;
1334     clevel2->trigty = -1.;
1335     clevel2->wartrig = 0.;
1336     clevel2->good = 0;
1337     clevel2->nstrip = 0.;
1338     clevel2->qtot = 0.;
1339     // clevel2->impx = 0.;
1340     // clevel2->impy = 0.;
1341     clevel2->tanx = 0.; // this is correct since it refers to the fortran structure
1342     clevel2->tany = 0.; // this is correct since it refers to the fortran structure
1343     clevel2->qmax = 0.;
1344     clevel2->nx22 = 0.;
1345     clevel2->qx22 = 0.;
1346     memset(clevel2->perr, 0, 4*sizeof(Int_t));
1347     memset(clevel2->swerr, 0, 4*sizeof(Int_t));
1348     memset(clevel2->crc, 0, 4*sizeof(Int_t));
1349     memset(clevel2->qq, 0, 4*sizeof(Int_t));
1350     memset(clevel2->varcfit, 0, 4*sizeof(Float_t));
1351     memset(clevel2->npcfit, 0, 4*sizeof(Int_t));
1352     memset(clevel2->planemax, 0, 2*sizeof(Int_t));
1353 mocchiut 1.6 memset(clevel2->selfdelay, 0, 4*7*sizeof(Int_t));
1354 mocchiut 1.1 memset(clevel2->fmode, 0, 2*sizeof(Int_t));
1355     memset(clevel2->cibar, 0, 2*22*sizeof(Int_t));
1356     memset(clevel2->cbar, 0, 2*22*sizeof(Float_t));
1357     }
1358    
1359     void CaloLevel0::ClearCalibVals(Int_t s){
1360     //
1361     for ( Int_t d=0 ; d<11 ;d++ ){
1362     Int_t pre = -1;
1363     for ( Int_t j=0; j<96 ;j++){
1364     if ( j%16 == 0 ) pre++;
1365     if ( s == 2 ){
1366     calped[0][2*d+1][j] = 0.;
1367     cstwerr[3] = 0.;
1368     cperror[3] = 0.;
1369     calgood[0][2*d+1][j] = 0.;
1370     calthr[0][2*d+1][pre] = 0.;
1371     calrms[0][2*d+1][j] = 0.;
1372     calbase[0][2*d+1][pre] = 0.;
1373     calvar[0][2*d+1][pre] = 0.;
1374     };
1375     if ( s == 3 ){
1376     calped[0][2*d][j] = 0.;
1377     cstwerr[1] = 0.;
1378     cperror[1] = 0.;
1379     calgood[0][2*d][j] = 0.;
1380     calthr[0][2*d][pre] = 0.;
1381     calrms[0][2*d][j] = 0.;
1382     calbase[0][2*d][pre] = 0.;
1383     calvar[0][2*d][pre] = 0.;
1384     };
1385     if ( s == 0 ){
1386     calped[1][2*d][j] = 0.;
1387     cstwerr[0] = 0.;
1388     cperror[0] = 0.;
1389     calgood[1][2*d][j] = 0.;
1390     calthr[1][2*d][pre] = 0.;
1391     calrms[1][2*d][j] = 0.;
1392     calbase[1][2*d][pre] = 0.;
1393     calvar[1][2*d][pre] = 0.;
1394     };
1395     if ( s == 1 ){
1396     calped[1][2*d+1][j] = 0.;
1397     cstwerr[2] = 0.;
1398     cperror[2] = 0.;
1399     calgood[1][2*d+1][j] = 0.;
1400     calthr[1][2*d+1][pre] = 0.;
1401     calrms[1][2*d+1][j] = 0.;
1402     calbase[1][2*d+1][pre] = 0.;
1403     calvar[1][2*d+1][pre] = 0.;
1404     };
1405     };
1406     };
1407     return;
1408     }
1409    
1410 mocchiut 1.7 Int_t CaloLevel0::Update(GL_TABLES *glt, UInt_t atime, Int_t s){
1411 mocchiut 1.1 //
1412 mocchiut 1.7 const TString host = glt->CGetHost();
1413     const TString user = glt->CGetUser();
1414     const TString psw = glt->CGetPsw();
1415     TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
1416     if ( !dbc->IsConnected() ) throw -116;
1417 mocchiut 1.8 stringstream myquery;
1418     myquery.str("");
1419     myquery << "SET time_zone='+0:00'";
1420     dbc->Query(myquery.str().c_str());
1421 mocchiut 1.1 Int_t sgnl = 0;
1422     //
1423     GL_CALO_CALIB *glcalo = new GL_CALO_CALIB();
1424     //
1425     sgnl = 0;
1426     //
1427     idcalib[s] = 0;
1428     fromtime[s] = 0;
1429     totime[s] = 0;
1430     calibno[s] = 0;
1431     ClearCalibVals(s);
1432     //
1433     UInt_t uptime = 0;
1434     //
1435     sgnl = glcalo->Query_GL_CALO_CALIB(atime,uptime,s,dbc);
1436     if ( sgnl < 0 ){
1437     if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
1438     return(sgnl);
1439     };
1440     //
1441     idcalib[s] = glcalo->ID_ROOT_L0;
1442     fromtime[s] = glcalo->FROM_TIME;
1443     if ( glcalo->TO_TIME < atime ){ // calibration is corrupted and we are using the one that preceed the good one
1444     totime[s] = uptime;
1445     } else {
1446     totime[s] = glcalo->TO_TIME;
1447     };
1448     // totime[s] = glcalo->TO_TIME;
1449     calibno[s] = glcalo->EV_ROOT;
1450     //
1451     if ( totime[s] == 0 ){
1452     if ( verbose ) printf(" CALORIMETER - WARNING: data with no associated calibration\n");
1453     ClearCalibVals(s);
1454     sgnl = 100;
1455     };
1456     //
1457     // determine path and name and entry of the calibration file
1458     //
1459     GL_ROOT *glroot = new GL_ROOT();
1460     if ( verbose ) printf("\n");
1461     if ( verbose ) printf(" ** SECTION %i **\n",s);
1462     //
1463     sgnl = glroot->Query_GL_ROOT(idcalib[s],dbc);
1464     if ( sgnl < 0 ){
1465     if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
1466     return(sgnl);
1467     };
1468     //
1469     stringstream name;
1470     name.str("");
1471     name << glroot->PATH.Data() << "/";
1472     name << glroot->NAME.Data();
1473     //
1474     fcalname[s] = (TString)name.str().c_str();
1475     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]);
1476     //
1477     sgnl = LoadCalib(s);
1478     //
1479     if ( sgnl != 0 ) return(sgnl);
1480     delete glcalo;
1481     delete glroot;
1482     //
1483     return(0);
1484     //
1485     }
1486    
1487     Int_t CaloLevel0::LoadCalib(Int_t s){
1488     //
1489     ifstream myfile;
1490     myfile.open(fcalname[s].Data());
1491     if ( !myfile ){
1492     return(-107);
1493     };
1494     myfile.close();
1495     //
1496     TFile *File = new TFile(fcalname[s].Data());
1497     if ( !File ) return(-108);
1498     TTree *tr = (TTree*)File->Get("CalibCalPed");
1499     if ( !tr ) return(-109);
1500     //
1501     TBranch *calo = tr->GetBranch("CalibCalPed");
1502     //
1503     pamela::CalibCalPedEvent *ce = 0;
1504     tr->SetBranchAddress("CalibCalPed", &ce);
1505     //
1506     Long64_t ncalibs = calo->GetEntries();
1507     //
1508     if ( !ncalibs ) return(-110);
1509     //
1510     calo->GetEntry(calibno[s]);
1511     //
1512     if (ce->cstwerr[s] != 0 && ce->cperror[s] == 0 ) {
1513     for ( Int_t d=0 ; d<11 ;d++ ){
1514     Int_t pre = -1;
1515     for ( Int_t j=0; j<96 ;j++){
1516     if ( j%16 == 0 ) pre++;
1517     if ( s == 2 ){
1518     calped[0][2*d+1][j] = ce->calped[3][d][j];
1519     cstwerr[3] = ce->cstwerr[3];
1520     cperror[3] = ce->cperror[3];
1521     calgood[0][2*d+1][j] = ce->calgood[3][d][j];
1522     calthr[0][2*d+1][pre] = ce->calthr[3][d][pre];
1523     calrms[0][2*d+1][j] = ce->calrms[3][d][j];
1524     calbase[0][2*d+1][pre] = ce->calbase[3][d][pre];
1525     calvar[0][2*d+1][pre] = ce->calvar[3][d][pre];
1526     };
1527     if ( s == 3 ){
1528     calped[0][2*d][j] = ce->calped[1][d][j];
1529     cstwerr[1] = ce->cstwerr[1];
1530     cperror[1] = ce->cperror[1];
1531     calgood[0][2*d][j] = ce->calgood[1][d][j];
1532     calthr[0][2*d][pre] = ce->calthr[1][d][pre];
1533     calrms[0][2*d][j] = ce->calrms[1][d][j];
1534     calbase[0][2*d][pre] = ce->calbase[1][d][pre];
1535     calvar[0][2*d][pre] = ce->calvar[1][d][pre];
1536     };
1537     if ( s == 0 ){
1538     calped[1][2*d][j] = ce->calped[0][d][j];
1539     cstwerr[0] = ce->cstwerr[0];
1540     cperror[0] = ce->cperror[0];
1541     calgood[1][2*d][j] = ce->calgood[0][d][j];
1542     calthr[1][2*d][pre] = ce->calthr[0][d][pre];
1543     calrms[1][2*d][j] = ce->calrms[0][d][j];
1544     calbase[1][2*d][pre] = ce->calbase[0][d][pre];
1545     calvar[1][2*d][pre] = ce->calvar[0][d][pre];
1546     };
1547     if ( s == 1 ){
1548     calped[1][2*d+1][j] = ce->calped[2][d][j];
1549     cstwerr[2] = ce->cstwerr[2];
1550     cperror[2] = ce->cperror[2];
1551     calgood[1][2*d+1][j] = ce->calgood[2][d][j];
1552     calthr[1][2*d+1][pre] = ce->calthr[2][d][pre];
1553     calrms[1][2*d+1][j] = ce->calrms[2][d][j];
1554     calbase[1][2*d+1][pre] = ce->calbase[2][d][pre];
1555     calvar[1][2*d+1][pre] = ce->calvar[2][d][pre];
1556     };
1557     };
1558     };
1559     } else {
1560     if ( verbose ) printf(" CALORIMETER - ERROR: problems finding a good calibration in this file! \n\n ");
1561     return(-111);
1562     };
1563     File->Close();
1564     return(0);
1565     }

  ViewVC Help
Powered by ViewVC 1.1.23