/[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.22 - (hide annotations) (download)
Thu Jun 19 20:05:44 2008 UTC (16 years, 5 months ago) by mocchiut
Branch: MAIN
Changes since 1.21: +891 -392 lines
New calorimeter calibration

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 mocchiut 1.22 memset(ctsicor, 0, 2*22*9*sizeof(Float_t));
71 mocchiut 1.9 memset(ctneigcor, 0, 2*22*6*sizeof(Float_t));
72 mocchiut 1.1 calopar1 = true;
73     calopar2 = true;
74     calopar3 = true;
75 mocchiut 1.22 calopar4 = true;
76     calopar5 = true;
77 mocchiut 1.1 crosst = true;
78     ftcalopar1 = 0;
79     ttcalopar1 = 0;
80     ftcalopar2 = 0;
81     ttcalopar2 = 0;
82     ftcalopar3 = 0;
83     ttcalopar3 = 0;
84 mocchiut 1.22 ftcalopar4 = 0;
85     ttcalopar4 = 0;
86     ftcalopar5 = 0;
87     ttcalopar5 = 0;
88 mocchiut 1.1 }
89    
90     void CaloLevel0::SetCrossTalk(Bool_t ct){
91     crosst = ct;
92 mocchiut 1.12 }
93 mocchiut 1.1
94 mocchiut 1.9 void CaloLevel0::SetCrossTalkType(Bool_t ct){
95     ctground = ct;
96 mocchiut 1.12 }
97 mocchiut 1.9
98 mocchiut 1.22 void CaloLevel0::SetCrossTalkType(Int_t ct){
99     if ( ct == 0 ) ctground = true;
100     if ( ct == 1 ){
101     ctground = false;
102     noselfct = false;
103     };
104     if ( ct == 2 ){
105     ctground = false;
106     noselfct = true;
107     };
108     }
109    
110 mocchiut 1.1 void CaloLevel0::SetVerbose(Bool_t ct){
111     verbose = ct;
112 mocchiut 1.12 }
113 mocchiut 1.1
114     /**
115     * Initialize CaloLevel0 object
116     **/
117 mocchiut 1.21 void CaloLevel0::ProcessingInit(TSQLServer *dbc, UInt_t hs, Int_t &sgnl, TTree *l0tree, Bool_t isdeb, Bool_t isverb){
118     if ( !dbc->IsConnected() ) throw -116;
119     this->InitDo(dbc,hs,sgnl,l0tree,isdeb,isverb);
120     }
121    
122     /**
123     * Initialize CaloLevel0 object
124     **/
125 mocchiut 1.7 void CaloLevel0::ProcessingInit(GL_TABLES *glt, UInt_t hs, Int_t &sgnl, TTree *l0tree, Bool_t isdeb, Bool_t isverb){
126     //
127     const TString host = glt->CGetHost();
128     const TString user = glt->CGetUser();
129     const TString psw = glt->CGetPsw();
130     TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
131     if ( !dbc->IsConnected() ) throw -116;
132 mocchiut 1.21 this->InitDo(dbc,hs,sgnl,l0tree,isdeb,isverb);
133     dbc->Close();
134     delete dbc;
135     }
136    
137    
138     void CaloLevel0::InitDo(TSQLServer *dbc, UInt_t hs, Int_t &sgnl, TTree *l0tree, Bool_t isdeb, Bool_t isverb){
139 mocchiut 1.8 stringstream myquery;
140     myquery.str("");
141     myquery << "SET time_zone='+0:00'";
142     dbc->Query(myquery.str().c_str());
143 mocchiut 1.1 //
144     debug = isdeb;
145     verbose = isverb;
146     //
147     l0tr=(TTree*)l0tree;
148     de = new pamela::calorimeter::CalorimeterEvent();
149     l0calo = (TBranch*)l0tr->GetBranch("Calorimeter");
150     l0tr->SetBranchAddress("Calorimeter", &de);
151     //
152     trkseqno = 0;
153     ClearStructs();
154     //
155     GL_CALO_CALIB *glcalo = new GL_CALO_CALIB();
156     //
157     sgnl = 0;
158     UInt_t uptime = 0;
159     //
160     for (Int_t s = 0; s < 4; s++){
161     idcalib[s] = 0;
162     fromtime[s] = 0;
163     totime[s] = 0;
164     calibno[s] = 0;
165     ClearCalibVals(s);
166     //
167     sgnl = glcalo->Query_GL_CALO_CALIB(hs,uptime,s,dbc);
168     if ( sgnl < 0 ){
169     if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
170     return;
171     };
172     //
173     idcalib[s] = glcalo->ID_ROOT_L0;
174     fromtime[s] = glcalo->FROM_TIME;
175     if ( glcalo->TO_TIME < hs ){ // calibration is corrupted and we are using the one that preceed the good one
176     totime[s] = uptime;
177     } else {
178     totime[s] = glcalo->TO_TIME;
179     };
180     calibno[s] = glcalo->EV_ROOT;
181     //
182     if ( totime[s] == 0 ){
183     if ( verbose ) printf(" CALORIMETER - WARNING: data with no associated calibration\n");
184     ClearCalibVals(s);
185     sgnl = 100;
186     };
187     };
188     //
189     // determine path and name and entry of the calibration file
190     //
191     GL_ROOT *glroot = new GL_ROOT();
192     if ( verbose ) printf("\n");
193     for (Int_t s = 0; s < 4; s++){
194     if ( verbose ) printf(" ** SECTION %i **\n",s);
195     if ( totime[s] > 0 ){
196     //
197     sgnl = glroot->Query_GL_ROOT(idcalib[s],dbc);
198     if ( sgnl < 0 ){
199     if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
200     return;
201     };
202     //
203     stringstream name;
204     name.str("");
205     name << glroot->PATH.Data() << "/";
206     name << glroot->NAME.Data();
207     //
208     fcalname[s] = (TString)name.str().c_str();
209     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]);
210     } else {
211     if ( verbose ) printf(" - runheader at time %u. NO CALIBRATION INCLUDE THE RUNHEADER! ",hs);
212     };
213     sgnl = LoadCalib(s);
214     if ( sgnl ) break;
215     };
216     //
217     delete glcalo;
218     delete glroot;
219     //
220     return;
221     //
222     }
223    
224 mocchiut 1.7 Int_t CaloLevel0::ChkCalib(GL_TABLES *glt, UInt_t atime){
225 mocchiut 1.1 Int_t sgnl = 0;
226     for ( Int_t s = 0; s < 4; s++){
227     if ( atime > totime[s] ){
228 mocchiut 1.7 sgnl = Update(glt,atime,s);
229 mocchiut 1.1 if ( sgnl < 0 ) return(sgnl);
230     };
231     };
232     return(sgnl);
233     }
234    
235 mocchiut 1.21 Int_t CaloLevel0::ChkParam(TSQLServer *dbc, UInt_t runheader, Bool_t mechal){
236     Int_t sig = this->ChkParamDo(dbc,runheader,mechal);
237     return(sig);
238     }
239    
240 mocchiut 1.7 Int_t CaloLevel0::ChkParam(GL_TABLES *glt, UInt_t runheader, Bool_t mechal){
241     const TString host = glt->CGetHost();
242     const TString user = glt->CGetUser();
243     const TString psw = glt->CGetPsw();
244     TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
245     if ( !dbc->IsConnected() ) throw -116;
246 mocchiut 1.8 stringstream myquery;
247     myquery.str("");
248     myquery << "SET time_zone='+0:00'";
249     dbc->Query(myquery.str().c_str());
250 mocchiut 1.7 //
251 mocchiut 1.21 Int_t sig = this->ChkParamDo(dbc,runheader,mechal);
252     dbc->Close();
253     delete dbc;
254     return(sig);
255     }
256    
257     Int_t CaloLevel0::ChkParamDo(TSQLServer *dbc, UInt_t runheader, Bool_t mechal){
258     //
259 mocchiut 1.1 stringstream calfile;
260     stringstream bmfile;
261     stringstream aligfile;
262     Int_t error = 0;
263     FILE *f = 0;
264     ifstream badfile;
265     GL_PARAM *glparam = new GL_PARAM();
266     //
267     if ( calopar1 || ( ttcalopar1 != 0 && ttcalopar1 < runheader ) ){
268 mocchiut 1.2 //
269     if ( debug ) printf(" calopar1 %i ftcalopar1 %u ttcalopar1 %u runheader %u \n",calopar1,ftcalopar1,ttcalopar1,runheader);
270     //
271 mocchiut 1.22 if ( calopar1 ){
272     //
273     // determine where I can find calorimeter ADC to MIP conversion file
274     //
275     if ( verbose ) printf(" Querying DB for calorimeter parameters files...\n");
276     //
277     error = 0;
278     error = glparam->Query_GL_PARAM(runheader,101,dbc);
279     if ( error < 0 ) return(error);
280     //
281     calfile.str("");
282     calfile << glparam->PATH.Data() << "/";
283     calfile << glparam->NAME.Data();
284     //
285     if ( verbose ) printf("\n Using ADC to MIP conversion file: \n %s \n",calfile.str().c_str());
286     f = fopen(calfile.str().c_str(),"rb");
287     if ( !f ){
288     if ( verbose ) printf(" CALORIMETER - ERROR: no ADC to MIP file!\n");
289     return(-105);
290     };
291     //
292     for (Int_t m = 0; m < 2 ; m++ ){
293     for (Int_t k = 0; k < 22; k++ ){
294     for (Int_t l = 0; l < 96; l++ ){
295     fread(&mip[m][k][l],sizeof(mip[m][k][l]),1,f);
296     if ( debug ) printf(" %f \n",mip[m][k][l]);
297     };
298     };
299     };
300     fclose(f);
301     };
302     //
303 mocchiut 1.1 calopar1 = false;
304     //
305 mocchiut 1.22 // flight extra corrections:
306 mocchiut 1.1 //
307 mocchiut 1.22 if ( verbose ) printf(" Querying DB for calorimeter flight ADC to MIP files...\n");
308 mocchiut 1.1 //
309     error = 0;
310 mocchiut 1.22 error = glparam->Query_GL_PARAM(runheader,110,dbc);
311 mocchiut 1.1 if ( error < 0 ) return(error);
312     //
313     calfile.str("");
314     calfile << glparam->PATH.Data() << "/";
315     calfile << glparam->NAME.Data();
316     ftcalopar1 = glparam->FROM_TIME;
317     ttcalopar1 = glparam->TO_TIME;
318     //
319 mocchiut 1.22 if ( verbose ) printf("\n Using ADC to MIP special conversion file: \n %s \n",calfile.str().c_str());
320     ifstream spfile;
321     spfile.open(calfile.str().c_str());
322     if ( !spfile ){
323     if ( verbose ) printf(" CALORIMETER - ERROR: no special calibration file!\n");
324     return(-123);
325 mocchiut 1.1 };
326 mocchiut 1.22 //
327     Int_t vview = 0;
328     Int_t vplane = 0;
329     Int_t vstrip = 0;
330     Float_t vval = 0.;
331     while ( spfile >> vview && spfile >> vplane && spfile >> vstrip && spfile >> vval){
332     if ( debug ) printf(" Setting ADC to MIP conversion factor: view %i plane %i strip %i mip %f \n",vview,vplane,vstrip,vval);
333     mip[vview][vplane][vstrip] = vval;
334     };
335 mocchiut 1.1 //
336     };
337     //
338 mocchiut 1.22 //
339 mocchiut 1.1 if ( calopar2 || ( ttcalopar2 != 0 && ttcalopar2 < runheader ) ){
340 mocchiut 1.21 //
341 mocchiut 1.2 if ( debug ) printf(" calopar2 %i ftcalopar2 %u ttcalopar2 %u runheader %u \n",calopar2,ftcalopar2,ttcalopar2,runheader);
342 mocchiut 1.1 calopar2 = false;
343     //
344     // determine where I can find calorimeter alignment file
345     //
346     //
347     error = 0;
348     error = glparam->Query_GL_PARAM(runheader,102,dbc);
349     if ( error < 0 ) return(error);
350     //
351     aligfile.str("");
352     aligfile << glparam->PATH.Data() << "/";
353     aligfile << glparam->NAME.Data();
354     ftcalopar2 = glparam->FROM_TIME;
355     ttcalopar2 = glparam->TO_TIME;
356     //
357 mocchiut 1.2 if ( verbose ) printf("\n Using parameter file: \n %s \n",aligfile.str().c_str());
358 mocchiut 1.1 f = fopen(aligfile.str().c_str(),"rb");
359     if ( !f ){
360 mocchiut 1.2 if ( verbose ) printf(" CALORIMETER - ERROR: no parameter file!\n");
361 mocchiut 1.1 return(-106);
362     };
363     //
364 mocchiut 1.2 if ( !mechal ){
365     //
366     fread(&clevel1->xalig,sizeof(clevel1->xalig),1,f);
367     if ( debug ) printf(" xalig = %f \n",clevel1->xalig);
368     fread(&clevel1->yalig,sizeof(clevel1->yalig),1,f);
369     if ( debug ) printf(" yalig = %f \n",clevel1->yalig);
370     fread(&clevel1->zalig,sizeof(clevel1->zalig),1,f);
371     if ( debug ) printf(" zalig = %f \n",clevel1->zalig);
372     } else {
373     if ( verbose ) printf("\n Using MECHANICAL alignement parameters \n");
374     //
375     CaloStrip cs = CaloStrip();
376     cs.UseMechanicalAlig();
377     clevel1->xalig = cs.GetXalig();
378     if ( debug ) printf(" xalig = %f \n",clevel1->xalig);
379     clevel1->yalig = cs.GetYalig();
380     if ( debug ) printf(" yalig = %f \n",clevel1->yalig);
381     clevel1->zalig = cs.GetZalig();
382     if ( debug ) printf(" zalig = %f \n",clevel1->zalig);
383     //
384     Float_t tmp = 0;
385     fread(&tmp,sizeof(clevel1->xalig),1,f);
386     fread(&tmp,sizeof(clevel1->yalig),1,f);
387     fread(&tmp,sizeof(clevel1->zalig),1,f);
388 mocchiut 1.11 // clevel1->zalig = -265.82;
389 mocchiut 1.2 //
390     };
391 mocchiut 1.1 fread(&clevel1->emin,sizeof(clevel1->emin),1,f);
392     if ( debug ) printf(" signal threshold = %f \n",clevel1->emin);
393     //
394     fclose(f);
395     };
396     //
397     // Load offline bad strip mask
398     //
399     if ( calopar3 || ( ttcalopar3 != 0 && ttcalopar3 < runheader ) ){
400 mocchiut 1.2 if ( debug ) printf(" calopar3 %i ftcalopar3 %u ttcalopar3 %u runheader %u \n",calopar3,ftcalopar3,ttcalopar3,runheader);
401 mocchiut 1.1 calopar3 = false;
402     //
403     // determine where I can find calorimeter alignment file
404     //
405     //
406     error = 0;
407     error = glparam->Query_GL_PARAM(runheader,103,dbc);
408     if ( error < 0 ) return(error);
409     //
410     bmfile.str("");
411     bmfile << glparam->PATH.Data() << "/";
412     bmfile << glparam->NAME.Data();
413     ftcalopar3 = glparam->FROM_TIME;
414     ttcalopar3 = glparam->TO_TIME;
415     //
416     if ( verbose ) printf("\n Using bad strip offline mask file: \n %s \n\n",bmfile.str().c_str());
417     badfile.open(bmfile.str().c_str());
418     if ( !badfile ){
419     if ( verbose ) printf(" CALORIMETER - ERROR: no bad strip offline mask file!\n");
420     return(-115);
421     };
422     //
423     Bool_t isdone = false;
424     Int_t bad = 0;
425     Int_t view = 1;
426     Int_t strip = 0;
427     Int_t plane = 21;
428     while ( !isdone ) {
429     badfile >> bad;
430     obadmask[view][plane][strip] = bad;
431     if ( debug && bad ) printf(" SETTING view %i plane %i strip %i BAD = %i \n",view,plane,strip,bad);
432     strip++;
433     if ( strip > 95 ){
434     strip = 0;
435     plane--;
436     if ( plane < 0 ){
437     plane = 21;
438     view--;
439     };
440     if ( view < 0 ) isdone = true;
441     };
442     };
443     //
444     badfile.close();
445     };
446     //
447 mocchiut 1.22 // calopar4
448     //
449     if ( calopar4 || ( ttcalopar4 != 0 && ttcalopar4 < runheader ) ){
450     //
451     if ( debug ) printf(" calopar4 %i ftcalopar4 %u ttcalopar4 %u runheader %u \n",calopar4,ftcalopar4,ttcalopar4,runheader);
452     //
453     calopar4 = false;
454     //
455     // flight extra corrections:
456     //
457     if ( verbose ) printf(" Querying DB for calorimeter max rms file...\n");
458     //
459     error = 0;
460     error = glparam->Query_GL_PARAM(runheader,109,dbc);
461     if ( error < 0 ) return(error);
462     //
463     calfile.str("");
464     calfile << glparam->PATH.Data() << "/";
465     calfile << glparam->NAME.Data();
466     ftcalopar4 = glparam->FROM_TIME;
467     ttcalopar4 = glparam->TO_TIME;
468     //
469     if ( verbose ) printf("\n Using calorimeter max rms file: \n %s \n",calfile.str().c_str());
470     ifstream spfile;
471     spfile.open(calfile.str().c_str());
472     if ( !spfile ){
473     if ( verbose ) printf(" CALORIMETER - ERROR: no max rms file!\n");
474     return(-124);
475     };
476     //
477     Int_t vview = 0;
478     Int_t vplane = 0;
479     Int_t vval = 0;
480     for (Int_t l=0; l<2; l++){
481     for (Int_t m=0; m<22; m++){
482     maxrms[l][m] = 26;
483     };
484     };
485     while ( spfile >> vview && spfile >> vplane && spfile >> vval){
486     if ( debug ) printf(" Setting view %i plane %i max rms %i \n",vview,vplane,vval);
487     maxrms[vview][vplane] = vval;
488     };
489     spfile.close();
490     //
491     };
492     //
493     // calopar5
494     //
495     if ( calopar5 || ( ttcalopar5 != 0 && ttcalopar5 < runheader ) ){
496     //
497     if ( debug ) printf(" calopar5 %i ftcalopar5 %u ttcalopar5 %u runheader %u \n",calopar5,ftcalopar5,ttcalopar5,runheader);
498     //
499     calopar5 = false;
500     //
501     // flight extra corrections:
502     //
503     if ( verbose ) printf(" Querying DB for calorimeter noise to signal threshold file...\n");
504     //
505     error = 0;
506     error = glparam->Query_GL_PARAM(runheader,111,dbc);
507     if ( error < 0 ) return(error);
508     //
509     calfile.str("");
510     calfile << glparam->PATH.Data() << "/";
511     calfile << glparam->NAME.Data();
512     ftcalopar5 = glparam->FROM_TIME;
513     ttcalopar5 = glparam->TO_TIME;
514     //
515     if ( verbose ) printf("\n Using calorimeter noise to signal threshold file: \n %s \n",calfile.str().c_str());
516     ifstream spfile;
517     spfile.open(calfile.str().c_str());
518     if ( !spfile ){
519     if ( verbose ) printf(" CALORIMETER - ERROR: no noise to signal threshold file!\n");
520     return(-125);
521     };
522     //
523     Int_t vview = 0;
524     Int_t vplane = 0;
525     Int_t vstrip = 0;
526     Float_t vval = 0.;
527     for (Int_t l=0; l<2; l++){
528     for (Int_t m=0; m<22; m++){
529     for (Int_t n=0; n<96; n++){
530     memin[l][m][n] = 0.7;
531     };
532     };
533     };
534     while ( spfile >> vview && spfile >> vplane && spfile >> vstrip && spfile >> vval){
535     if ( vstrip == -1 ){
536     for (Int_t ll=0; ll<96; ll++){
537     if ( debug ) printf(" Setting view %i plane %i strip %i noise to signal ratio %f \n",vview,vplane,ll,vval);
538     memin[vview][vplane][ll] = vval;
539     };
540     } else {
541     if ( debug ) printf(" Setting view %i plane %i strip %i noise to signal ratio %f \n",vview,vplane,vstrip,vval);
542     memin[vview][vplane][vstrip] = vval;
543     };
544     };
545     spfile.close();
546     //
547     };
548     //
549     //
550 mocchiut 1.1 delete glparam;
551     //
552     return(0);
553     }
554    
555 mocchiut 1.22 Int_t CaloLevel0::CalcCrossTalkCorr(TSQLServer *dbc, UInt_t runheader, Bool_t ctusetable){
556     Int_t sig = CalcCrossTalkCorrDo(dbc,runheader,ctusetable);
557     return(sig);
558     };
559    
560 mocchiut 1.21 Int_t CaloLevel0::CalcCrossTalkCorr(TSQLServer *dbc, UInt_t runheader){
561 mocchiut 1.22 Int_t sig = CalcCrossTalkCorrDo(dbc,runheader,true);
562 mocchiut 1.21 return(sig);
563     }
564    
565 mocchiut 1.22 Int_t CaloLevel0::CalcCrossTalkCorr(GL_TABLES *glt, UInt_t runheader, Bool_t usetable){
566     const TString host = glt->CGetHost();
567     const TString user = glt->CGetUser();
568     const TString psw = glt->CGetPsw();
569     TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
570     if ( !dbc->IsConnected() ) throw -116;
571     stringstream myquery;
572     myquery.str("");
573     myquery << "SET time_zone='+0:00'";
574     dbc->Query(myquery.str().c_str());
575     //
576     Int_t sig = CalcCrossTalkCorrDo(dbc,runheader,usetable);
577     dbc->Close();
578     delete dbc;
579     //
580     return(sig);
581     //
582     };
583    
584 mocchiut 1.9 Int_t CaloLevel0::CalcCrossTalkCorr(GL_TABLES *glt, UInt_t runheader){
585     const TString host = glt->CGetHost();
586     const TString user = glt->CGetUser();
587     const TString psw = glt->CGetPsw();
588     TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
589     if ( !dbc->IsConnected() ) throw -116;
590     stringstream myquery;
591     myquery.str("");
592     myquery << "SET time_zone='+0:00'";
593     dbc->Query(myquery.str().c_str());
594     //
595 mocchiut 1.22 Int_t sig = CalcCrossTalkCorrDo(dbc,runheader,true);
596 mocchiut 1.21 dbc->Close();
597     delete dbc;
598     //
599     return(sig);
600     //
601     }
602    
603 mocchiut 1.22 Int_t CaloLevel0::CalcCrossTalkCorrDo(TSQLServer *dbc, UInt_t runheader, Bool_t usetable){
604 mocchiut 1.21 //
605     if ( ctground ) return(0);
606     //
607 mocchiut 1.9 Int_t error = 0;
608 mocchiut 1.22 GL_PARAM *glparam = new GL_PARAM();
609 mocchiut 1.9 //
610     // determine where I can find file with offline bad pulser mask
611     //
612 mocchiut 1.22 stringstream bmfile;
613 mocchiut 1.9 error = 0;
614     error = glparam->Query_GL_PARAM(runheader,105,dbc);
615     if ( error < 0 ) return(error);
616     //
617     bmfile.str("");
618     bmfile << glparam->PATH.Data() << "/";
619     bmfile << glparam->NAME.Data();
620     //
621 mocchiut 1.22 ifstream badfile;
622 mocchiut 1.9 if ( verbose ) printf("\n Using bad pulser offline mask file: \n %s \n\n",bmfile.str().c_str());
623     badfile.open(bmfile.str().c_str());
624     if ( !badfile ){
625     if ( verbose ) printf(" CALORIMETER - ERROR: no bad pulser offline mask file!\n");
626     return(-115);
627     };
628     //
629     Bool_t isdone = false;
630     Int_t bad = 0;
631     Int_t view = 1;
632     Int_t pre = 0;
633     Int_t plane = 21;
634     while ( !isdone ) {
635     badfile >> bad;
636     obadpulsemask[view][plane][pre] = bad;
637     if ( debug && bad ) printf(" SETTING view %i plane %i pre %i BAD = %i \n",view,plane,pre,bad);
638     pre++;
639     if ( pre > 5 ){
640     pre = 0;
641     plane--;
642     if ( plane < 0 ){
643     plane = 21;
644     view--;
645     };
646     if ( view < 0 ) isdone = true;
647     };
648     };
649     //
650     badfile.close();
651 mocchiut 1.22 if ( !usetable ){
652 mocchiut 1.9 //
653 mocchiut 1.22 // Let's start with cross-talk correction calculation
654 mocchiut 1.9 //
655 mocchiut 1.22 GL_CALOPULSE_CALIB *glp = new GL_CALOPULSE_CALIB();
656     Float_t adcp[2][22][96];
657     Float_t adcpcal[2][22][96];
658     memset(adcp , 0, 2*22*96*sizeof(Float_t));
659     memset(adcpcal , 0, 2*22*96*sizeof(Float_t));
660 mocchiut 1.9 //
661 mocchiut 1.22 UInt_t pampli = 0;
662     for (Int_t s=0; s<4; s++){
663     //
664     // Save into matrix adcp the values of the highest pulse calibration (pulse amplitude = 2)
665     //
666     pampli = 2;
667     error = 0;
668     error = glp->Query_GL_CALOPULSE_CALIB(runheader,s,pampli,dbc);
669     if ( error < 0 ){
670     if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
671     return(error);
672     };
673     //
674     UInt_t idcalib = glp->ID_ROOT_L0;
675     UInt_t fromtime = glp->FROM_TIME;
676     UInt_t calibno = glp->EV_ROOT;
677     //
678     // determine path and name and entry of the calibration file
679     //
680     GL_ROOT *glroot = new GL_ROOT();
681     if ( verbose ) printf("\n");
682     if ( verbose ) printf(" ** SECTION %i **\n",s);
683     //
684     error = 0;
685     error = glroot->Query_GL_ROOT(idcalib,dbc);
686     if ( error < 0 ){
687     if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
688     return(error);
689     };
690     //
691     stringstream name;
692     name.str("");
693     name << glroot->PATH.Data() << "/";
694     name << glroot->NAME.Data();
695     //
696     TString fcalname = (TString)name.str().c_str();
697     ifstream myfile;
698     myfile.open(fcalname.Data());
699     if ( !myfile ){
700     return(-107);
701     };
702     myfile.close();
703     //
704     TFile *File = new TFile(fcalname.Data());
705     if ( !File ) return(-108);
706     TTree *tr = (TTree*)File->Get("CalibCalPulse2");
707     if ( !tr ) return(-119);
708     //
709     TBranch *calo = tr->GetBranch("CalibCalPulse2");
710     //
711     pamela::CalibCalPulse2Event *ce = 0;
712     tr->SetBranchAddress("CalibCalPulse2", &ce);
713     //
714     Long64_t ncalibs = calo->GetEntries();
715     //
716     if ( !ncalibs ) return(-110);
717     //
718     calo->GetEntry(calibno);
719     if ( verbose ) printf(" PULSE2 using entry %u from file %s",calibno,fcalname.Data());
720     //
721     // retrieve calibration table
722     //
723     if ( ce->pstwerr[s] && ce->pperror[s] == 0 && ce->unpackError == 0 ){
724     for ( Int_t d=0 ; d<11 ;d++ ){
725     for ( Int_t j=0; j<96 ;j++){
726     if ( s == 2 ){
727     adcp[0][2*d+1][j] = ce->calpuls[3][d][j];
728     };
729     if ( s == 3 ){
730     adcp[0][2*d][j] = ce->calpuls[1][d][j];
731     };
732     if ( s == 0 ){
733     adcp[1][2*d][j] = ce->calpuls[0][d][j];
734     };
735     if ( s == 1 ){
736     adcp[1][2*d+1][j] = ce->calpuls[2][d][j];
737     };
738     };
739     };
740     } else {
741     if ( verbose ) printf(" CALORIMETER - ERROR: problems finding a good calibration in this file! \n\n ");
742     return(-111);
743     };
744     //
745     File->Close();
746     delete glroot;
747     //
748     // Save into matrix adcpcal the calibrated values of the pulse calibration (subtraction of pulse amplitude = 0 relative to the pulse2 calibration used)
749     //
750     pampli = 0;
751     error = 0;
752     error = glp->Query_GL_CALOPULSE_CALIB(fromtime,s,pampli,dbc);
753     if ( error < 0 ){
754     if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
755     return(error);
756     };
757     //
758     idcalib = glp->ID_ROOT_L0;
759     calibno = glp->EV_ROOT;
760     //
761     // determine path and name and entry of the calibration file
762     //
763     glroot = new GL_ROOT();
764     if ( verbose ) printf("\n");
765     if ( verbose ) printf(" ** SECTION %i **\n",s);
766     //
767     error = 0;
768     error = glroot->Query_GL_ROOT(idcalib,dbc);
769     if ( error < 0 ){
770     if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
771     return(error);
772     };
773     //
774     name.str("");
775     name << glroot->PATH.Data() << "/";
776     name << glroot->NAME.Data();
777     //
778     fcalname = (TString)name.str().c_str();
779     myfile.open(fcalname.Data());
780     if ( !myfile ){
781     return(-107);
782     };
783     myfile.close();
784     //
785     TFile *File1 = new TFile(fcalname.Data());
786     if ( !File1 ) return(-108);
787     TTree *tr1 = (TTree*)File1->Get("CalibCalPulse1");
788     if ( !tr1 ) return(-120);
789     //
790     TBranch *calo1 = tr1->GetBranch("CalibCalPulse1");
791     //
792     pamela::CalibCalPulse1Event *ce1 = 0;
793     tr1->SetBranchAddress("CalibCalPulse1", &ce1);
794     //
795     ncalibs = calo1->GetEntries();
796     //
797     if ( !ncalibs ) return(-110);
798     //
799     calo1->GetEntry(calibno);
800     if ( verbose ) printf(" PULSE1 using entry %u from file %s",calibno,fcalname.Data());
801     //
802     // retrieve calibration table
803     //
804     if ( ce1->pstwerr[s] && ce1->pperror[s] == 0 && ce1->unpackError == 0 ){
805     for ( Int_t d=0 ; d<11 ;d++ ){
806     for ( Int_t j=0; j<96 ;j++){
807     if ( s == 2 ){
808     adcpcal[0][2*d+1][j] = adcp[0][2*d+1][j] - ce1->calpuls[3][d][j];
809     };
810     if ( s == 3 ){
811     adcpcal[0][2*d][j] = adcp[0][2*d][j] - ce1->calpuls[1][d][j];
812     };
813     if ( s == 0 ){
814     adcpcal[1][2*d][j] = adcp[1][2*d][j] - ce1->calpuls[0][d][j];
815     };
816     if ( s == 1 ){
817     adcpcal[1][2*d+1][j] = adcp[1][2*d+1][j] - ce1->calpuls[2][d][j];
818     };
819     };
820     };
821     } else {
822     if ( verbose ) printf(" CALORIMETER - ERROR: problems finding a good calibration in this file! \n\n ");
823     return(-111);
824     };
825     //
826     File1->Close();
827     //
828     delete glroot;
829     //
830     };// loop on the four sections
831 mocchiut 1.9 //
832     //
833 mocchiut 1.22 delete glp;
834 mocchiut 1.9 //
835 mocchiut 1.22 // Ok, now we can try to calculate the cross-talk correction for each pre-amplifier
836 mocchiut 1.9 //
837 mocchiut 1.22 for ( Int_t v=0; v<2; v++){
838     if ( debug ) printf(" \n\n NEW VIEW \n");
839     for ( Int_t p=0; p<22; p++){
840     for ( Int_t npre=0; npre<6; npre++){
841     ctprecor[v][p][npre] = 1000.;
842     ctneigcor[v][p][npre] = 1000.;
843     Int_t str0=npre*16;
844     Int_t str16= -1 + (1+npre)*16;
845     //
846     UInt_t neigc = 0;
847     UInt_t farc = 0;
848     UInt_t pulsc = 0;
849     Float_t sigpulsed = 0.;
850     Float_t neigbase = 0.;
851     Float_t farbase = 0.;
852     //
853     // 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
854     // moreover count the number of strips in each case
855     //
856     for (Int_t s=str0; s<=str16; s++){
857     if ( adcpcal[v][p][s] > 10000.){
858     sigpulsed = adcpcal[v][p][s];
859     pulsc++;
860     if ( s > str0 ){
861     neigbase += adcpcal[v][p][s-1];
862     neigc++;
863     farbase -= adcpcal[v][p][s-1];
864     farc--;
865     };
866     if ( s < str16 ){
867     neigbase += adcpcal[v][p][s+1];
868     neigc++;
869     farbase -= adcpcal[v][p][s+1];
870     farc--;
871     };
872     } else {
873     farc++;
874     farbase += adcpcal[v][p][s];
875     };
876 mocchiut 1.9 };
877 mocchiut 1.22 //
878     // Now calculate the corrections
879     //
880     Float_t avefarbase = 0.;
881     if ( farc ) avefarbase = farbase/(Float_t)farc;
882     Float_t aveneigbase = 0.;
883     if ( neigc ) aveneigbase = neigbase/(Float_t)neigc;
884     //
885     if ( pulsc == 1 && farc && neigc ){
886     ctprecor[v][p][npre] = -avefarbase/(sigpulsed+fabs(avefarbase));
887     ctneigcor[v][p][npre] = fabs(aveneigbase-avefarbase)/(sigpulsed+fabs(avefarbase));
888     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]);
889     } else {
890     //
891     // did not find the pulsed strip or more than one pulsed strip found!
892     //
893     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);
894     //
895 mocchiut 1.9 };
896     };
897 mocchiut 1.22 if ( debug ) printf(" \n ==================== \n");
898 mocchiut 1.9 };
899     };
900 mocchiut 1.22 } else {
901 mocchiut 1.9 //
902 mocchiut 1.22 // use pre-amply table
903     //
904 mocchiut 1.9 //
905 mocchiut 1.22 // determine where I can find file with offline neighbour correction table
906 mocchiut 1.9 //
907 mocchiut 1.22 stringstream bmfile2;
908 mocchiut 1.9 error = 0;
909 mocchiut 1.22 error = glparam->Query_GL_PARAM(runheader,106,dbc);
910     if ( error < 0 ) return(error);
911     //
912     bmfile2.str("");
913     bmfile2 << glparam->PATH.Data() << "/";
914     bmfile2 << glparam->NAME.Data();
915     //
916     ifstream badfile2;
917     if ( verbose ) printf("\n Using pre-amply neighbour crosstalk table file: \n %s \n\n",bmfile2.str().c_str());
918     badfile2.open(bmfile2.str().c_str());
919     if ( !badfile2 ){
920     if ( verbose ) printf(" CALORIMETER - ERROR: no pre-amply neighbour crosstalk table file!\n");
921     return(-121);
922 mocchiut 1.9 };
923 mocchiut 1.22 //
924     Int_t vview = 0;
925     Int_t vplane = 0;
926     Int_t vpre = 0;
927     Float_t vcorr = 0.;
928     while ( badfile2 >> vview && badfile2 >> vplane && badfile2 >> vpre && badfile2 >> vcorr){
929     if ( debug ) printf(" Pre-amply neighbour correction: view %i plane %i pre %i correction %f \n",vview,vplane,vpre,vcorr);
930     ctneigcor[vview][vplane][vpre] = vcorr;
931     };
932 mocchiut 1.9 //
933 mocchiut 1.22 // determine where I can find file with offline SECOND neighbour correction table
934 mocchiut 1.9 //
935 mocchiut 1.22 stringstream bmfile3;
936 mocchiut 1.9 error = 0;
937 mocchiut 1.22 error = glparam->Query_GL_PARAM(runheader,107,dbc);
938     if ( error < 0 ) return(error);
939 mocchiut 1.9 //
940 mocchiut 1.22 bmfile3.str("");
941     bmfile3 << glparam->PATH.Data() << "/";
942     bmfile3 << glparam->NAME.Data();
943     //
944     ifstream badfile3;
945     if ( verbose ) printf("\n Using pre-amply second neighbour crosstalk table file: \n %s \n\n",bmfile3.str().c_str());
946     badfile3.open(bmfile3.str().c_str());
947     if ( !badfile3 ){
948     if ( verbose ) printf(" CALORIMETER - ERROR: no pre-amply second neighbour crosstalk table file!\n");
949     return(-122);
950 mocchiut 1.9 };
951 mocchiut 1.22 //
952     Int_t pview = 0;
953     Int_t pplane = 0;
954     Int_t ppre = 0;
955     Float_t pcorr = 0.;
956     while ( badfile3 >> pview && badfile3 >> pplane && badfile3 >> ppre && badfile3 >> pcorr){
957     if ( debug ) printf(" Pre-amply second neighbour correction: view %i plane %i pre %i correction %f \n",pview,pplane,ppre,-pcorr);
958     ctprecor[pview][pplane][ppre] = -pcorr; // data are saved as negatives in the file
959     };
960 mocchiut 1.9 //
961 mocchiut 1.22 // determine where to find the file containing the Silicon crosstalk correction table
962 mocchiut 1.9 //
963 mocchiut 1.22 stringstream bmfile4;
964     error = 0;
965     error = glparam->Query_GL_PARAM(runheader,108,dbc);
966     if ( error < 0 ) return(error);
967 mocchiut 1.9 //
968 mocchiut 1.22 bmfile4.str("");
969     bmfile4 << glparam->PATH.Data() << "/";
970     bmfile4 << glparam->NAME.Data();
971     //
972     ifstream badfile4;
973     if ( verbose ) printf("\n Using Silicon crosstalk table file: \n %s \n\n",bmfile4.str().c_str());
974     badfile4.open(bmfile4.str().c_str());
975     if ( !badfile4 ){
976     if ( verbose ) printf(" CALORIMETER - ERROR: no Silicon crosstalk table file!\n");
977     return(-125);
978 mocchiut 1.9 };
979 mocchiut 1.22 //
980     Int_t spview = 0;
981     Int_t spplane = 0;
982     Int_t psil = 0;
983     Float_t spcorr = 0.;
984     memset(ctsicor, 0, 2*22*9*sizeof(Float_t));
985     while ( badfile4 >> spview && badfile4 >> spplane && badfile4 >> psil && badfile4 >> spcorr){
986     if ( debug ) printf(" Silicon correction: view %i plane %i silicon %i correction %f \n",spview,spplane,psil,-spcorr);
987     ctsicor[spview][spplane][psil] = -spcorr; // data are saved as negatives in the file
988     };
989 mocchiut 1.9 //
990 mocchiut 1.22 };
991 mocchiut 1.9 //
992 mocchiut 1.22 delete glparam;
993 mocchiut 1.9 //
994     // Check the calculated corrections
995     //
996     Int_t opre=0;
997     Int_t ppre=0;
998     Bool_t found = false;
999     for ( Int_t v=0; v<2; v++){
1000     for ( Int_t p=0; p<22; p++){
1001     for ( Int_t npre=0; npre<6; npre++){
1002     if ( ctprecor[v][p][npre] == 1000. || ctneigcor[v][p][npre] == 1000. || obadpulsemask[v][p][npre] != 0 ){
1003     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]);
1004     if ( npre%2 ){
1005     opre = npre-1;
1006     } else {
1007     opre = npre+1;
1008     };
1009     if ( ctprecor[v][p][opre] == 1000. || ctneigcor[v][p][opre] == 1000. || obadpulsemask[v][p][opre] != 0 ){
1010     ppre=0;
1011     found = false;
1012     while ( ppre < 6 ){
1013     if ( ctprecor[v][p][ppre] != 1000. && ctneigcor[v][p][ppre] != 1000. && !obadpulsemask[v][p][ppre] ){
1014     found = true;
1015     ctprecor[v][p][npre] = ctprecor[v][p][ppre];
1016     ctneigcor[v][p][npre] = ctneigcor[v][p][ppre];
1017     break;
1018     };
1019     ppre++;
1020     };
1021     if ( !found ){
1022     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);
1023     ctprecor[v][p][npre] = 0.002;
1024     ctneigcor[v][p][npre] = 0.002;
1025     };
1026     } else {
1027     ctprecor[v][p][npre] = ctprecor[v][p][opre];
1028     ctneigcor[v][p][npre] = ctneigcor[v][p][opre];
1029     };
1030     if ( debug ) printf(" AFTER: pre-correction: %f neighbour strips correction %f \n",ctprecor[v][p][npre],ctneigcor[v][p][npre]);
1031     };
1032     };
1033     };
1034 mocchiut 1.22 };
1035 mocchiut 1.9 //
1036     return(0);
1037 mocchiut 1.12 }
1038 mocchiut 1.1
1039 mocchiut 1.16 void CaloLevel0::FindBaseCompress(Int_t l, Int_t m, Int_t pre){
1040 mocchiut 1.17 Int_t n = 0;
1041     Float_t q = 0;
1042     this->FindBaseCompress(l,m,pre,n,q);
1043     }
1044    
1045     void CaloLevel0::FindBaseCompress(Int_t l, Int_t m, Int_t pre, Int_t &nst, Float_t &qp){
1046 mocchiut 1.16 for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
1047     dexy[l][m][e] = dexyc[l][m][e];
1048     };
1049 mocchiut 1.17 this->FindBaseRaw(l,m,pre,nst,qp);
1050 mocchiut 1.16 }
1051    
1052 mocchiut 1.1 void CaloLevel0::FindBaseRaw(Int_t l, Int_t m, Int_t pre){
1053 mocchiut 1.17 Int_t n = 0;
1054     Float_t q = 0;
1055     this->FindBaseRaw(l,m,pre,n,q);
1056     }
1057    
1058     void CaloLevel0::FindBaseRaw(Int_t l, Int_t m, Int_t pre, Int_t &nst, Float_t &qp){
1059     //
1060     Float_t minstrip = 100000.;
1061     Float_t rms = 0.;
1062 mocchiut 1.21 Int_t process = 0;
1063     Int_t onlmask[16];
1064     memset(onlmask, 0, 16*sizeof(Int_t));
1065     //
1066     while ( process < 2 ){
1067     //
1068     minstrip = 100000.;
1069     rms = 0.;
1070     base[l][m][pre] = 0.;
1071     qp = 0.;
1072     //
1073     Int_t spos = -1;
1074     Int_t ee = 0;
1075 mocchiut 1.1 for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
1076 mocchiut 1.21 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. && onlmask[ee] == 0 ) {
1077     minstrip = dexy[l][m][e]-calped[l][m][e];
1078     rms = calthr[l][m][pre];
1079     spos = ee;
1080 mocchiut 1.17 };
1081 mocchiut 1.21 ee++;
1082     qp += (dexy[l][m][e]-calped[l][m][e]-sbase[l][m][e]);
1083 mocchiut 1.1 };
1084 mocchiut 1.17 //
1085 mocchiut 1.21 if ( debug && l==0 ){
1086     printf("\n BASELINE CALCULATION for view %i pl %i pre %i: \n => minstrip %f rms %f \n => qp = %f \n",l,m,pre,minstrip,rms,qp);
1087 mocchiut 1.16 };
1088 mocchiut 1.21 if ( minstrip != 100000. ) {
1089     Float_t strip6s = 0.;
1090     for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
1091     if ( (dexy[l][m][e]-calped[l][m][e]) >= minstrip && (dexy[l][m][e]-calped[l][m][e]) <= (minstrip+rms) ) {
1092     strip6s += 1.;
1093     base[l][m][pre] += (dexy[l][m][e] - calped[l][m][e]);
1094     };
1095     //
1096     // compression
1097     //
1098     // if ( abs((int)(dexy[l][m][e]-calped[l][m][e])) <= (minstrip+rms) ) {
1099     // dexyc[l][m][e] = 0.;
1100     // } else {
1101     dexyc[l][m][e] = dexy[l][m][e];
1102     // };
1103 mocchiut 1.17 };
1104     //
1105 mocchiut 1.21 if ( strip6s == 1. && process < 1 ){
1106     onlmask[spos] = 1;
1107     process++;
1108 mocchiut 1.22 if ( debug ) printf(" Warning, only one strip to calculate baseline: minstrip %f rms %f spos %i l %i m %i pre %i \n",minstrip,rms,spos,l,m,pre);
1109 mocchiut 1.21 continue;
1110     };
1111     process += 2;
1112     nst = (Int_t)strip6s;
1113 mocchiut 1.17 //
1114 mocchiut 1.22 if ( debug ){
1115 mocchiut 1.21 printf(" strip6s %f \n",strip6s);
1116     };
1117     // if ( strip6s >= 9. ){
1118     if ( (strip6s >= 2. && process == 2) || (strip6s >= 9. && process > 2) ){
1119 mocchiut 1.22 //if ( (strip6s >= 4. && process == 2) || (strip6s >= 9. && process > 2) ){
1120 mocchiut 1.21 Double_t arro = base[l][m][pre]/strip6s;
1121     Float_t deci = 1000.*((float)arro - float(int(arro)));
1122     if ( deci < 500. ) {
1123     arro = double(int(arro));
1124     } else {
1125     arro = 1. + double(int(arro));
1126     };
1127     base[l][m][pre] = arro;
1128     //
1129     // if too few strips were used to determine the baseline check if it is comparable with the previous event, if not mark it as bad
1130     //
1131 mocchiut 1.22 if ( debug && process > 2 ) printf(" AGH low strip value was discarded process %i strip6s %f minstrip %f rms %f spos %i\n",process,strip6s,minstrip,rms,spos);
1132     if ( debug ) printf(" Calculated baseline: base %f sbase-0.02*qp %f \n",base[l][m][pre],(-qp*0.02+sbase[l][m][pre]));
1133 mocchiut 1.21 //
1134 mocchiut 1.22 if ( strip6s < 4 && base[l][m][pre] > (-0.015*qp+sbase[l][m][pre]) && sbase[l][m][pre] > 0. ){
1135     if ( debug ) printf(" Suspicious calculated baseline: base %f sbase-0.02*qp %f strip6s %i \n",base[l][m][pre],(-qp*0.02+sbase[l][m][pre]),(Int_t)strip6s);
1136 mocchiut 1.21 base[l][m][pre] = 31000.;
1137     for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
1138     dexyc[l][m][e] = dexy[l][m][e];
1139     };
1140     };
1141     } else {
1142     base[l][m][pre] = 31000.;
1143     for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
1144     dexyc[l][m][e] = dexy[l][m][e];
1145     };
1146 mocchiut 1.17 };
1147 mocchiut 1.1 } else {
1148 mocchiut 1.21 process += 2;
1149 mocchiut 1.17 base[l][m][pre] = 31000.;
1150     for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
1151     dexyc[l][m][e] = dexy[l][m][e];
1152     };
1153 mocchiut 1.1 };
1154 mocchiut 1.17 };
1155 mocchiut 1.1 }
1156    
1157     Int_t CaloLevel0::Calibrate(Int_t ei){
1158     //
1159     // get entry ei
1160     //
1161     l0calo->GetEntry(ei);
1162     //
1163     // if it was not a selftrigger event, could it ever been a selftrigger event? if so trigty = 3.
1164     //
1165 mocchiut 1.14 clevel2->nsatstrip = 0.;
1166 mocchiut 1.1 Int_t val = 0;
1167     Int_t del = 1100;
1168 mocchiut 1.6 for (Int_t sec = 0; sec < 4; sec++){
1169     for (Int_t dsec = 0; dsec < 7; dsec++){
1170     val = (Int_t)de->calselftrig[sec][dsec];
1171     del = delay(val);
1172     clevel2->selfdelay[sec][dsec] = del;
1173     };
1174     };
1175     val = 0;
1176     del = 1100;
1177 mocchiut 1.1 if ( clevel2->trigty != 2. ){
1178     Bool_t ck = false;
1179     for (Int_t sec = 0; sec < 4; sec++){
1180     val = (Int_t)de->calselftrig[sec][6];
1181     del = delay(val);
1182     if ( del < 1100 ){
1183     clevel2->wartrig = 0.;
1184     clevel2->trigty = 3.;
1185     ck = true;
1186     break;
1187     };
1188     };
1189     if ( !ck ) clevel2->wartrig = 100.;
1190     } else {
1191     Bool_t ck = false;
1192     for (Int_t sec = 0; sec < 4; sec++){
1193     val = (Int_t)de->calselftrig[sec][6];
1194     del = delay(val);
1195     if ( del < 1100 ){
1196     clevel2->wartrig = 0.;
1197     ck = true;
1198     };
1199     };
1200     if ( !ck ) clevel2->wartrig = 100.;
1201     };
1202     //
1203     Int_t se = 5;
1204     Int_t done = 0;
1205     Int_t pre = -1;
1206     Bool_t isCOMP = false;
1207     Bool_t isFULL = false;
1208     Bool_t isRAW = false;
1209     Float_t ener;
1210     Int_t doneb = 0;
1211     Int_t donec = 0;
1212 mocchiut 1.22 Int_t ck[2][22][6];
1213     memset(ck, 0, 2*22*6*sizeof(Int_t));
1214 mocchiut 1.1 Int_t ipre = 0;
1215 mocchiut 1.19 // Int_t ip[3] = {0};
1216     Int_t ip[3] = {0,0,0};
1217 mocchiut 1.22 Int_t ipp = 0;
1218 mocchiut 1.1 Float_t base0, base1, base2;
1219     base0 = 0.;
1220     base1 = 0.;
1221     base2 = 0.;
1222 mocchiut 1.22 Float_t qpre[2][22][6];
1223     memset(qpre, 0, 2*22*6*sizeof(Float_t));
1224 mocchiut 1.1 Float_t ene[96];
1225     Int_t chdone[4] = {0,0,0,0};
1226     Int_t pe = 0;
1227     //
1228     Float_t ener0 = 0.;
1229     Float_t cbase0 = 0.;
1230     Bool_t pproblem = false;
1231 mocchiut 1.22 Bool_t negbase = false;
1232 mocchiut 1.1 //
1233     Float_t tim = 0.;
1234     Int_t plo = 0;
1235     Int_t fbi = 0;
1236     Int_t cle = 0;
1237     //
1238     // run over views and planes
1239     //
1240     for (Int_t l = 0; l < 2; l++){
1241     for (Int_t m = 0; m < 22; m++){
1242     //
1243 mocchiut 1.2 // determine the section number
1244 mocchiut 1.1 //
1245 mocchiut 1.22 negbase = false;
1246 mocchiut 1.1 se = 5;
1247 mocchiut 1.2 if (l == 0 && m%2 == 0) se = 3;
1248 mocchiut 1.1 if (l == 0 && m%2 != 0) se = 2;
1249 mocchiut 1.2 if (l == 1 && m%2 != 0) se = 1;
1250     if (l == 1 && m%2 == 0) se = 0;
1251 mocchiut 1.1 //
1252     // determine what kind of event we are going to analyze
1253     //
1254     isCOMP = false;
1255     isFULL = false;
1256     isRAW = false;
1257     if ( de->stwerr[se] & (1 << 16) ) isCOMP = true;
1258     if ( de->stwerr[se] & (1 << 17) ) isFULL = true;
1259     if ( de->stwerr[se] & (1 << 3) ) isRAW = true;
1260     if ( !chdone[se] ){
1261     //
1262     // check for any error in the event
1263     //
1264     clevel2->crc[se] = 0;
1265     if ( de->perror[se] == 132 ){
1266     clevel2->crc[se] = 1;
1267     pe++;
1268     };
1269     clevel2->perr[se] = 0;
1270     if ( de->perror[se] != 0 ){
1271     clevel2->perr[se] = 1;
1272     pe++;
1273     };
1274     clevel2->swerr[se] = 0;
1275     for (Int_t j = 0; j < 7 ; j++){
1276     if ( (j != 3) && (de->stwerr[se] & (1 << j)) ){
1277     clevel2->swerr[se] = 1;
1278     pe++;
1279     };
1280     };
1281     chdone[se] = 1;
1282     };
1283     if ( clevel2->crc[se] == 0 && (clevel1->good2 == 1 || clevel2->trigty >= 2) ){
1284     pre = -1;
1285     //
1286     for (Int_t nn = 0; nn < 96; nn++){
1287 mocchiut 1.22 // ene[nn] = 0.;
1288 mocchiut 1.1 dexy[l][m][nn] = de->dexy[l][m][nn] ;
1289     dexyc[l][m][nn] = de->dexyc[l][m][nn] ;
1290     };
1291     //
1292     // run over preamplifiers
1293     //
1294     pre = -1;
1295     cbase0 = 0.;
1296 mocchiut 1.17 Int_t nstt[2];
1297     Float_t rqp[2];
1298 mocchiut 1.1 for (Int_t i = 0; i < 3; i++){
1299 mocchiut 1.21 nstt[0] = 1000;
1300     nstt[1] = 1000;
1301 mocchiut 1.17 rqp[0] = 0.;
1302     rqp[1] = 0.;
1303 mocchiut 1.1 for (Int_t j = 0; j < 2; j++){
1304     pre = j + i*2;
1305     //
1306     // baseline check and calculation
1307     //
1308 mocchiut 1.16 if ( !isRAW ){
1309 mocchiut 1.18 //
1310     // if it is a compress event with fully transmitted pre try to calculate the baseline
1311     //
1312 mocchiut 1.16 if ( de->base[l][m][pre] != 0. && de->base[l][m][pre]<31000. ) {
1313     base[l][m][pre] = de->base[l][m][pre] ;
1314     } else {
1315 mocchiut 1.17 FindBaseCompress(l,m,pre,nstt[j],rqp[j]);
1316 mocchiut 1.16 };
1317 mocchiut 1.1 cbase0 += base[l][m][pre];
1318     } else {
1319     //
1320 mocchiut 1.17 // if it is a raw event calculate the baseline.
1321 mocchiut 1.1 //
1322 mocchiut 1.17 FindBaseRaw(l,m,pre,nstt[j],rqp[j]);
1323 mocchiut 1.1 cbase0 += base[l][m][pre];
1324     };
1325     };
1326 mocchiut 1.17 //
1327     // if we are able to calculate the baseline with more than 3 strips on one pre and not in the other one choose the pre with more calculated strips
1328     //
1329 mocchiut 1.21 if ( nstt[0] < 4 && nstt[1] >= 4 && nstt[0] != 1000 && nstt[1] != 1000 ) base[l][m][pre-1] = 31000.;
1330     if ( nstt[0] >= 4 && nstt[1] < 4 && nstt[0] != 1000 && nstt[1] != 1000 ) base[l][m][pre] = 31000.;
1331 mocchiut 1.22 // //
1332     // // if we are NOT able to calculate the baseline with more than 3 strips on both pres take the baseline (if any) of the one which has less energy
1333     // //
1334     // if ( nstt[0] < 4 && nstt[1] < 4 ){
1335     // if ( rqp[0] >= rqp[1] ) base[l][m][pre-1] = 31000.;
1336     // if ( rqp[0] < rqp[1] ) base[l][m][pre] = 31000.;
1337     // };
1338 mocchiut 1.1 };
1339     //
1340     // run over strips
1341     //
1342     pre = -1;
1343     ener0 = 0.;
1344     for (Int_t i = 0 ; i < 3 ; i++){
1345     ip[i] = 0;
1346     for (Int_t n = i*32 ; n < (i+1)*32 ; n++){
1347     if (n%16 == 0) {
1348     done = 0;
1349     doneb = 0;
1350     donec = 0;
1351     pre++;
1352 mocchiut 1.22 ck[l][m][pre] = 0;
1353     qpre[l][m][pre] = 0.;
1354 mocchiut 1.1 };
1355     //
1356     // baseline check and calculation
1357     //
1358     // no suitable new baseline, use old ones!
1359     //
1360     if ( !done ){
1361     if ( (base[l][m][pre] == 31000. || base[l][m][pre] == 0.) ){
1362 mocchiut 1.22 ck[l][m][pre] = 1;
1363 mocchiut 1.1 if (pre%2 == 0) {
1364     ip[i] = pre + 1;
1365     } else {
1366     ip[i] = pre - 1;
1367     };
1368 mocchiut 1.3 if ( (base[l][m][ip[i]] == 31000. || base[l][m][ip[i]] == 0. || !crosst ) ){
1369 mocchiut 1.1 //
1370 mocchiut 1.22 ck[l][m][pre] = 2;
1371 mocchiut 1.1 if ( sbase[l][m][pre] == 31000. || sbase[l][m][pre] == 0. ) {
1372 mocchiut 1.22 ck[l][m][pre] = 3;
1373 mocchiut 1.1 };
1374     };
1375     };
1376 mocchiut 1.21 done = 1;
1377 mocchiut 1.1 };
1378     //
1379     // CALIBRATION ALGORITHM
1380     //
1381     if ( !doneb ){
1382 mocchiut 1.22 if ( debug ) printf(" ck[l][m][pre] is %i \n",ck[l][m][pre]);
1383     switch (ck[l][m][pre]) {
1384 mocchiut 1.1 case 0:
1385     base0 = base[l][m][pre];
1386     base2 = calbase[l][m][pre];
1387 mocchiut 1.16 if ( debug ) printf(" base0 = base l%i m%i pre%i = %f base2 = calbase l m pre = %f \n",l,m,pre,base[l][m][pre],calbase[l][m][pre]);
1388 mocchiut 1.1 break;
1389     case 1:
1390     base0 = base[l][m][ip[i]];
1391     base2 = calbase[l][m][ip[i]];
1392 mocchiut 1.16 if ( debug ) printf(" base0 = base l%i m%i ip(i)%i = %f base2 = calbase l m ip(i) = %f \n",l,m,ip[i],base[l][m][ip[i]],calbase[l][m][ip[i]]);
1393 mocchiut 1.1 break;
1394     case 2:
1395     base0 = sbase[l][m][pre];
1396 mocchiut 1.3 base2 = calbase[l][m][pre];
1397 mocchiut 1.16 if ( debug ) printf(" base0 = sbase l%i m%i pre%i = %f base2 = calbase l m pre = %f \n",l,m,pre,sbase[l][m][pre],calbase[l][m][pre]);
1398 mocchiut 1.1 break;
1399     case 3:
1400     base0 = calbase[l][m][pre];
1401     base2 = calbase[l][m][pre];
1402 mocchiut 1.16 if ( debug ) printf(" base0 = calbase l%i m%i pre%i = %f base2 = calbase l m pre = %f \n",l,m,pre,calbase[l][m][pre],calbase[l][m][pre]);
1403 mocchiut 1.1 break;
1404     };
1405     base1 = calbase[l][m][pre];
1406     doneb = 1;
1407     };
1408     ener = dexyc[l][m][n];
1409     ener0 += ener;
1410     clevel1->estrip[n][m][l] = 0.;
1411 mocchiut 1.22 if ( de->base[l][m][pre] < 0 ) negbase = true;
1412 mocchiut 1.1 if ( base0>0 && base0 < 30000. ){
1413 mocchiut 1.22 //
1414     // save the baseline only if the energy release is "small"
1415     //
1416     if ( !donec && (base0 + base1 - base2) != 0. && (n+1)%16==0 ){
1417     if ( qpre[l][m][pre] < 200. ) sbase[l][m][pre] = base0 + base1 - base2;
1418 mocchiut 1.1 donec = 1;
1419     };
1420     if ( ener > 0. ){
1421     clevel1->estrip[n][m][l] = (ener - calped[l][m][n] - base0 - base1 + base2)/mip[l][m][n] ;
1422     //
1423     // 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)
1424     //
1425 mocchiut 1.22 if ( clevel1->estrip[n][m][l] > 0. ) qpre[l][m][pre] += clevel1->estrip[n][m][l];
1426 mocchiut 1.3 //
1427     //
1428 mocchiut 1.1 };
1429     };
1430     };
1431 mocchiut 1.22 };
1432     //
1433     // check if there were problems with 5.7 or glitches in the power supply
1434     //
1435     if ( ((ener0 == 0. && cbase0 == 0.) || negbase ) && !pproblem && clevel2->perr[se] == 0){
1436     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);
1437     pproblem = true;
1438     pe++;
1439     };
1440     //
1441     } else {
1442     for (Int_t nn = 0; nn < 96; nn++){
1443     clevel1->estrip[nn][m][l] = 0.;
1444     };
1445     };
1446     };
1447     };
1448     //
1449     // run over views and planes to apply crosstalk corrections
1450     //
1451     for (Int_t l = 0; l < 2; l++){
1452     for (Int_t m = 0; m < 22; m++){
1453     //
1454     // determine the section number
1455     //
1456     se = 5;
1457     if (l == 0 && m%2 == 0) se = 3;
1458     if (l == 0 && m%2 != 0) se = 2;
1459     if (l == 1 && m%2 != 0) se = 1;
1460     if (l == 1 && m%2 == 0) se = 0;
1461     //
1462     // check for any error in the event
1463     //
1464     if ( clevel2->crc[se] == 0 && (clevel1->good2 == 1 || clevel2->trigty >= 2) ){
1465     //
1466     // Cross-talk corrections
1467     //
1468     if ( crosst ){
1469     //
1470     // energy on silicon ladders
1471     //
1472     Float_t qsi[3];
1473     qsi[0] = qpre[l][m][0]+qpre[l][m][1];
1474     qsi[1] = qpre[l][m][2]+qpre[l][m][3];
1475     qsi[2] = qpre[l][m][4]+qpre[l][m][5];
1476     //
1477     for ( pre = 1; pre < 6; pre += 2 ){
1478     Int_t ladder = (pre - 1)/2;
1479     //
1480     // If the noselfct flag is set the strip doesn't suffer the self crosstalk due to electronics so we must subtract some energy
1481     //
1482     if ( noselfct ){
1483     for (Int_t j = ladder*32 ; j < (ladder+1)*32 ; j++){
1484     ipre = j/16 ;
1485     if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] -= clevel1->estrip[j][m][l] * ctprecor[l][m][ipre];
1486     };
1487     };
1488     //
1489     // Using the neighbour pre baseline
1490     //
1491     if (ck[l][m][pre] == 1 || ck[l][m][pre-1] == 1){
1492     //
1493     // pre-amplifier effect on baseline when using the neighbour pre (ck=1)
1494     //
1495     if (ck[l][m][pre] == 1){
1496 mocchiut 1.19 ipre = pre;
1497 mocchiut 1.22 ipp = pre - 1;
1498 mocchiut 1.1 } else {
1499 mocchiut 1.19 ipre = pre - 1;
1500 mocchiut 1.22 ipp = pre;
1501 mocchiut 1.1 };
1502 mocchiut 1.22 Int_t it = 0;
1503     Float_t nqpre = 0.;
1504     //
1505     if ( debug ) printf(" CK1 Limit for while: 0.07 \n");
1506 mocchiut 1.1 for (Int_t j = ipre*16 ; j < (ipre+1)*16 ; j++){
1507 mocchiut 1.9 if ( !ctground ){
1508 mocchiut 1.22 if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] += - qpre[l][m][ipp] * ctprecor[l][m][ipp];
1509 mocchiut 1.9 } else {
1510 mocchiut 1.22 if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] += - qpre[l][m][ipp] * 0.00478;
1511 mocchiut 1.9 };
1512 mocchiut 1.22 if ( clevel1->estrip[j][m][l] > 0. ) nqpre += clevel1->estrip[j][m][l] ;
1513 mocchiut 1.1 };
1514 mocchiut 1.22 qpre[l][m][ipre] = nqpre;
1515     nqpre = 0.;
1516     Float_t deltaqpre = qpre[l][m][ipre];
1517     //
1518     // these values are empirically determined, usually the routine converge due to deltaqsi and the latest applied correction is based on less than 1 mip
1519     //
1520     while ( it < 10 && deltaqpre > 0.07 ){
1521     nqpre = 0.;
1522     for (Int_t j = ipre*16 ; j < (ipre+1)*16 ; j++){
1523     if ( !ctground ){
1524     if ( debug ) printf(" CK1 pre correction: iteration %i deltaqpre %f ctprecor %f TOTAL CORRECTION %f \n",it,deltaqpre,ctprecor[l][m][ipre],deltaqpre * ctprecor[l][m][ipre]);
1525     if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] += deltaqpre * ctprecor[l][m][ipre];
1526     } else {
1527     if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] += deltaqpre * 0.00478;
1528     };
1529     if ( clevel1->estrip[j][m][l] > 0. ) nqpre += clevel1->estrip[j][m][l] ;
1530     };
1531     if ( ctground ) it = 100;
1532     it++;
1533     deltaqpre = nqpre - qpre[l][m][ipre];
1534     if ( debug ) printf(" CK1 BEFORE: qpre %f \n",qpre[l][m][ipre]);
1535     qpre[l][m][ipre] = nqpre;
1536     if ( debug ) printf(" CK1 AFTER: qpre %f \n",qpre[l][m][ipre]);
1537     };
1538     //
1539 mocchiut 1.1 };
1540 mocchiut 1.22 //
1541     // No baseline calculation due to high energy release
1542     //
1543     if (ck[l][m][pre] == 2 && ck[l][m][pre-1] == 2){
1544     //
1545     // y^
1546     // |
1547     // | 6 7 8
1548     // | 3 4 5
1549     // | 0 1 2
1550     // | --------------------------------------> x
1551     //
1552     Int_t si1 = 0;
1553     Int_t si2 = 0;
1554     Int_t si3 = 0;
1555     if ( l == 0 ){
1556     if ( ladder == 0 ){
1557     si1 = 0;
1558     si2 = 3;
1559     si3 = 6;
1560     };
1561     if ( ladder == 1 ){
1562     si1 = 1;
1563     si2 = 4;
1564     si3 = 7;
1565     };
1566     if ( ladder == 2 ){
1567     si1 = 2;
1568     si2 = 5;
1569     si3 = 8;
1570     };
1571     } else {
1572     if ( ladder == 0 ){
1573     si1 = 0;
1574     si2 = 1;
1575     si3 = 2;
1576     };
1577     if ( ladder == 1 ){
1578     si1 = 3;
1579     si2 = 4;
1580     si3 = 5;
1581     };
1582     if ( ladder == 2 ){
1583     si1 = 6;
1584     si2 = 7;
1585     si3 = 8;
1586     };
1587     };
1588     //
1589     // Find the energy distribution along the considered plane looking at the two sandwiching plane of the other view.
1590     //
1591     Float_t sied[3] = {0.,0.,0.};
1592     Int_t othv = !l;
1593     Int_t othpl1 = m - 1;
1594     Int_t othpl2 = m + 1;
1595     Float_t oprof[3] = {0.,0.,0.};
1596     for(Int_t s=0; s<3; s++){
1597     for(Int_t t=(s*32); t<32*(s + 1); t++){
1598     if ( othpl1 > -1 ) {
1599     oprof[s] += clevel1->estrip[othv][othpl1][t];
1600     };
1601     if ( othpl2 < 22 ) {
1602     oprof[s] += clevel1->estrip[othv][othpl2][t];
1603     };
1604     };
1605     };
1606     Float_t otote = fabs(oprof[0]) + fabs(oprof[1]) + fabs(oprof[2]);
1607     for(Int_t g=0; g<3; g++){
1608     if ( otote > 0. ){
1609     sied[g] = fabs(oprof[g])/otote;
1610 mocchiut 1.9 } else {
1611 mocchiut 1.22 sied[g] = 1./3.;
1612 mocchiut 1.9 };
1613 mocchiut 1.1 };
1614 mocchiut 1.22 //
1615     //
1616     //
1617     Int_t it = 0;
1618     Int_t jpre = 0;
1619     Float_t nqsi = 0.;
1620     Float_t snqsi = qsi[ladder];
1621     Float_t nqpre[2] = {0.,0.};
1622     Float_t deltaqsi = qsi[ladder];
1623     Float_t deltaqpre[2];
1624     deltaqpre[0] = qpre[l][m][pre-1];
1625     deltaqpre[1] = qpre[l][m][pre];
1626     //
1627     if ( debug ) printf(" Limit for while: 0.07 it < 10 \n");
1628     //
1629     // these values are empirically determined, usually the routine converge due to deltaqsi and the latest applied correction is based on less than 1 mip
1630     //
1631     while ( it < 10 && (deltaqsi > 0.07 || deltaqpre[0] > 0.07 || deltaqpre[1] > 0.07) ){
1632     nqsi = 0.;
1633     nqpre[0] = 0.;
1634     nqpre[1] = 0.;
1635     for (Int_t j = ladder*32 ; j < (ladder+1)*32 ; j++){
1636     ipre = 0;
1637     if ( j > (ladder*32)+15 ) ipre = 1;
1638     jpre = j/16 ;
1639     //
1640     // Silicon effect on the baseline when using the same pre previous baseline (ck = 2) + pre-amply effect
1641     //
1642     if ( !ctground ){
1643     if ( debug ) printf(" silicon correction: iteration %i deltaqsi[%i] %f ctsicor %f %f %f sied %f %f %f si %i %i %i TOTAL CORRECTION %f \n",it,ladder,deltaqsi,ctsicor[l][m][si1],ctsicor[l][m][si2],ctsicor[l][m][si3],sied[0],sied[1],sied[2],si1,si2,si3,deltaqsi * (ctsicor[l][m][si1] * sied[0] + ctsicor[l][m][si2] * sied[1] + ctsicor[l][m][si3] * sied[2]));
1644     if ( debug ) printf(" pre correction: iteration %i deltaqpre[0] %f deltaqpre[1] %f ctprecor %f TOTAL CORRECTION %f \n",it,deltaqpre[0],deltaqpre[1],ctprecor[l][m][jpre],deltaqpre[ipre] * ctprecor[l][m][jpre]);
1645     if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] += (deltaqsi * (ctsicor[l][m][si1] * sied[0] + ctsicor[l][m][si2] * sied[1] + ctsicor[l][m][si3] * sied[2])/mip[l][m][j]) + deltaqpre[ipre] * ctprecor[l][m][jpre];
1646     } else {
1647     if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] += 0. + qpre[l][m][jpre] * 0.00478; // no correction
1648     };
1649     if ( clevel1->estrip[j][m][l] > 0. ) nqsi += clevel1->estrip[j][m][l] ;
1650     if ( clevel1->estrip[j][m][l] > 0. ) nqpre[ipre] += clevel1->estrip[j][m][l] ;
1651     };
1652     if ( ctground ) it = 100;
1653     deltaqsi = nqsi-snqsi;
1654     snqsi = nqsi;
1655     it++;
1656     deltaqpre[0] = nqpre[0] - qpre[l][m][pre-1];
1657     deltaqpre[1] = nqpre[1] - qpre[l][m][pre];
1658     if ( debug ) printf(" BEFORE: qpre 0 %f qpre 1 %f \n",qpre[l][m][pre-1],qpre[l][m][pre]);
1659     qpre[l][m][pre-1] = nqpre[0];
1660     qpre[l][m][pre] = nqpre[1];
1661     if ( debug ) printf(" AFTER: qpre 0 %f qpre 1 %f \n",qpre[l][m][pre-1],qpre[l][m][pre]);
1662     };
1663     //
1664     //
1665     //
1666     // for (Int_t j = ladder*32 ; j < (ladder+1)*32 ; j++){
1667     // ipre = j/16 ;
1668     // //
1669     // // pre-amplifier effect on baseline when using the same pre previous event baseline (ck=2)
1670     // //
1671     // if ( !ctground ){
1672     // if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] += qpre[l][m][ipre] * ctprecor[l][m][ipre];
1673     // } else {
1674     // if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] += qpre[l][m][ipre] * 0.00478;
1675     // };
1676     // };
1677 mocchiut 1.1 };
1678     };
1679     };
1680 mocchiut 1.22 };
1681     //
1682     Int_t j4 = -4;
1683     Int_t jjj = -3;
1684     Int_t jj = -2;
1685     Int_t jjpre = -1;
1686     Int_t jjjpre = -1;
1687     memset(ene, 0, 96*sizeof(Float_t));
1688     for (Int_t j = 0 ; j < 100 ; j++){
1689     jj++;
1690     jjj++;
1691     j4++;
1692     if ( j < 96 ) ene[j] = clevel1->estrip[j][m][l];
1693     if ( crosst ){
1694     //
1695     // "Real" crosstalk effect on the neighbour strips respect to the one which have seen the energy deposit
1696     //
1697     if ( jj >= 0 && jj < 96 ){
1698     if ( !ctground ){
1699     if ( jj%16 == 0 ) jjpre++;
1700     if ( jj != 0 && jj != 32 && jj != 64 && ene[jj-1] != 0. ) ene[jj-1] += -clevel1->estrip[jj][m][l] * ctneigcor[l][m][jjpre];
1701     if ( jj != 31 && jj != 63 && jj != 95 && ene[jj+1] != 0. ) ene[jj+1] += -clevel1->estrip[jj][m][l] * ctneigcor[l][m][jjpre];
1702     } else {
1703     if ( jj != 0 && jj != 32 && jj != 64 && ene[jj-1] != 0. ) ene[jj-1] += -clevel1->estrip[jj][m][l] * 0.01581;
1704     if ( jj != 31 && jj != 63 && jj != 95 && ene[jj+1] != 0. ) ene[jj+1] += -clevel1->estrip[jj][m][l] * 0.01581;
1705 mocchiut 1.1 };
1706 mocchiut 1.22 };
1707     if ( jjj >= 0 && jjj < 96 ){
1708     if ( !ctground ){
1709     if ( jjj%16 == 0 ) jjjpre++;
1710     if ( jjj != 0 && jjj != 32 && jjj != 64 && clevel1->estrip[jjj-1][m][l] != 0. ) clevel1->estrip[jjj-1][m][l] += -ene[jjj] * ctneigcor[l][m][jjjpre];
1711     if ( jjj != 31 && jjj != 63 && jjj != 95 && clevel1->estrip[jjj+1][m][l] !=0. ) clevel1->estrip[jjj+1][m][l] += -ene[jjj] * ctneigcor[l][m][jjjpre];
1712     } else {
1713     if ( jjj != 0 && jjj != 32 && jjj != 64 && clevel1->estrip[jjj-1][m][l] != 0. ) clevel1->estrip[jjj-1][m][l] += -ene[jjj] * 0.01581;
1714     if ( jjj != 31 && jjj != 63 && jjj != 95 && clevel1->estrip[jjj+1][m][l] != 0. ) clevel1->estrip[jjj+1][m][l] += -ene[jjj] * 0.01581;
1715 mocchiut 1.1 };
1716     };
1717 mocchiut 1.22 };
1718     if ( j4 >= 0 && j4 < 96 ){
1719     //
1720     // CALOLEVEL1 CODING AND FILLING
1721     //
1722     //
1723     // NOTICE: THE FOLLOWING LINE EXCLUDE ALL STRIPS FOR WHICH THE RMS*4 IS GREATER THAN 26 !!! <=============== IMPORTANT! =================>
1724     //
1725     if ( obadmask[l][m][j4] == 1 || clevel1->estrip[j4][m][l] <= clevel1->emin || clevel1->estrip[j4][m][l] <= memin[l][m][j4] || calrms[l][m][j4] > maxrms[l][m] ){
1726     clevel1->estrip[j4][m][l] = 0.;
1727     };
1728     //
1729     // code and save the energy for each strip in svstrip
1730     //
1731     if ( clevel1->estrip[j4][m][l] > clevel1->emin ){
1732 mocchiut 1.1 //
1733 mocchiut 1.22 Float_t savel1 = clevel1->estrip[j4][m][l];
1734     // if ( dexyc[l][m][j4] == 32767. ){
1735     if ( dexyc[l][m][j4] > 32000. ){
1736     savel1 += 5000.;
1737     clevel2->nsatstrip += 1.;
1738 mocchiut 1.1 };
1739     //
1740 mocchiut 1.22 tim = 100000.;
1741     plo = m;
1742     fbi = 0;
1743     if ( savel1 > 0.99995 ){
1744     tim = 10000.;
1745 mocchiut 1.1 plo = m;
1746 mocchiut 1.22 fbi = 1;
1747     };
1748     if ( savel1 > 9.9995 ){
1749     tim = 1000.;
1750     plo = 22 + m;
1751     fbi = 1;
1752     };
1753     if ( savel1 > 99.995 ){
1754     tim = 100.;
1755     plo = 22 + m;
1756     fbi = 0;
1757     };
1758     if ( savel1 > 999.95 ){
1759     tim = 10.;
1760     plo = 44 + m;
1761     fbi = 0;
1762     };
1763     if ( savel1 > 9999.5 ){
1764     tim = 1.;
1765     plo = 66 + m;
1766 mocchiut 1.1 fbi = 0;
1767 mocchiut 1.22 };
1768     //
1769     cle = (Int_t)lroundf(tim*savel1);
1770     //
1771     if ( l == 0 ){
1772 mocchiut 1.1 //
1773 mocchiut 1.22 // +-PPSSmmmm.mmmm
1774 mocchiut 1.1 //
1775 mocchiut 1.22 svstrip[istrip] = fbi*1000000000 + plo*10000000 + j4*100000 + cle;
1776     } else {
1777     svstrip[istrip] = -(fbi*1000000000 + plo*10000000 + j4*100000 + cle);
1778 mocchiut 1.1 };
1779 mocchiut 1.22 //
1780     istrip++;
1781 mocchiut 1.1 };
1782     };
1783     };
1784 mocchiut 1.22 //
1785 mocchiut 1.1 };
1786     };
1787 mocchiut 1.22 //
1788     // store goodness flag
1789     //
1790 mocchiut 1.1 if ( !pe ){
1791     clevel2->good = 1;
1792     } else {
1793     clevel2->good = 0;
1794     };
1795 mocchiut 1.22 //
1796     // done
1797     //
1798 mocchiut 1.1 return(0);
1799     }
1800    
1801     void CaloLevel0::GetTrkVar(){
1802     calol2tr();
1803     }
1804    
1805     void CaloLevel0::FillTrkVar(CaloLevel2 *ca, Int_t nutrk){
1806     //
1807     CaloTrkVar *t_ca = new CaloTrkVar();
1808     //
1809     t_ca->trkseqno = trkseqno;
1810     t_ca->ncore = (Int_t)clevel2->ncore;
1811     t_ca->qcore = clevel2->qcore;
1812     t_ca->noint = (Int_t)clevel2->noint;
1813     t_ca->ncyl = (Int_t)clevel2->ncyl;
1814     t_ca->qcyl = clevel2->qcyl;
1815     t_ca->qtrack = clevel2->qtrack;
1816     t_ca->qtrackx = clevel2->qtrackx;
1817     t_ca->qtracky = clevel2->qtracky;
1818     t_ca->dxtrack = clevel2->dxtrack;
1819     t_ca->dytrack = clevel2->dytrack;
1820     t_ca->qlast = clevel2->qlast;
1821     t_ca->nlast = (Int_t)clevel2->nlast;
1822     t_ca->qpre = clevel2->qpre;
1823     t_ca->npre = (Int_t)clevel2->npre;
1824     t_ca->qpresh = clevel2->qpresh;
1825     t_ca->npresh = (Int_t)clevel2->npresh;
1826     t_ca->qtr = clevel2->qtr;
1827     t_ca->ntr = (Int_t)clevel2->ntr;
1828     t_ca->planetot = (Int_t)clevel2->planetot;
1829     t_ca->qmean = clevel2->qmean;
1830     t_ca->dX0l = clevel2->dX0l;
1831     t_ca->qlow = clevel2->qlow;
1832     t_ca->nlow = (Int_t)clevel2->nlow;
1833     //
1834     if ( trkseqno == -1 ){
1835     // ca->impx = clevel2->impx;
1836     // ca->impy = clevel2->impy;
1837     ca->tanx[1] = clevel2->tanx;
1838     ca->tany[1] = clevel2->tany;
1839     ca->elen = clevel2->elen;
1840     ca->selen = clevel2->selen;
1841     // memcpy(ca->cibar,clevel2->cibar,sizeof(clevel2->cibar));
1842     // memcpy(ca->cbar,clevel2->cbar,sizeof(clevel2->cbar));
1843     memcpy(t_ca->tibar,clevel2->cibar,sizeof(clevel2->cibar));
1844     memcpy(t_ca->tbar,clevel2->cbar,sizeof(clevel2->cbar));
1845     memcpy(ca->planemax,clevel2->planemax,sizeof(clevel2->planemax));
1846 mocchiut 1.6 memcpy(ca->selfdelay,clevel2->selfdelay,sizeof(clevel2->selfdelay));
1847 mocchiut 1.1 ca->varcfit[2] = clevel2->varcfit[0];
1848     ca->varcfit[3] = clevel2->varcfit[1];
1849     ca->npcfit[2] = clevel2->npcfit[0];
1850     ca->npcfit[3] = clevel2->npcfit[1];
1851     // memcpy(ca->varcfit,clevel2->varcfit,sizeof(clevel2->varcfit));
1852     // memcpy(ca->npcfit,clevel2->npcfit,sizeof(clevel2->npcfit));
1853     } else {
1854     memcpy(t_ca->tibar,clevel2->tibar,sizeof(clevel2->tibar));
1855     memcpy(t_ca->tbar,clevel2->tbar,sizeof(clevel2->tbar));
1856     };
1857     //
1858     //
1859     if ( !(ca->CaloTrk) ) ca->CaloTrk = new TClonesArray("CaloTrkVar",1); //ELENA
1860     TClonesArray &t = *ca->CaloTrk;
1861     new(t[nutrk]) CaloTrkVar(*t_ca);
1862     //
1863     delete t_ca;
1864     //
1865     ClearTrkVar();
1866     }
1867    
1868     void CaloLevel0::GetCommonVar(){
1869     calol2cm();
1870     }
1871    
1872     void CaloLevel0::FillCommonVar(CaloLevel1 *c1, CaloLevel2 *ca){
1873     //
1874     ca->good = clevel2->good;
1875     if ( clevel2->trigty == 2. ){
1876     ca->selftrigger = 1;
1877     } else {
1878     ca->selftrigger = 0;
1879     };
1880     //
1881     ca->selftrigger += (Int_t)clevel2->wartrig;
1882     //
1883     memcpy(ca->perr,clevel2->perr,sizeof(clevel2->perr));
1884     memcpy(ca->swerr,clevel2->swerr,sizeof(clevel2->swerr));
1885     memcpy(ca->crc,clevel2->crc,sizeof(clevel2->crc));
1886     ca->nstrip = (Int_t)clevel2->nstrip;
1887 mocchiut 1.13 ca->nsatstrip = (Int_t)clevel2->nsatstrip;
1888 mocchiut 1.1 ca->qtot = clevel2->qtot;
1889     // ca->impx = clevel2->impx;
1890     // ca->impy = clevel2->impy;
1891     ca->tanx[0] = clevel2->tanx;
1892     ca->tany[0] = clevel2->tany;
1893     ca->nx22 = (Int_t)clevel2->nx22;
1894     ca->qx22 = clevel2->qx22;
1895     ca->qmax = clevel2->qmax;
1896     ca->elen = clevel2->elen;
1897     ca->selen = clevel2->selen;
1898     memcpy(ca->qq,clevel2->qq,sizeof(clevel2->qq));
1899     memcpy(ca->planemax,clevel2->planemax,sizeof(clevel2->planemax));
1900 mocchiut 1.6 memcpy(ca->selfdelay,clevel2->selfdelay,sizeof(clevel2->selfdelay));
1901 mocchiut 1.1 ca->varcfit[0] = clevel2->varcfit[0];
1902     ca->varcfit[1] = clevel2->varcfit[1];
1903     ca->npcfit[0] = clevel2->npcfit[0];
1904     ca->npcfit[1] = clevel2->npcfit[1];
1905     ca->fitmode[0] = clevel2->fmode[0];
1906     ca->fitmode[1] = clevel2->fmode[1];
1907     // memcpy(ca->varcfit,clevel2->varcfit,sizeof(clevel2->varcfit));
1908     // memcpy(ca->npcfit,clevel2->npcfit,sizeof(clevel2->npcfit));
1909     memcpy(ca->cibar,clevel2->cibar,sizeof(clevel2->cibar));
1910     memcpy(ca->cbar,clevel2->cbar,sizeof(clevel2->cbar));
1911     //
1912     if ( c1 ){
1913     c1->istrip = istrip;
1914     c1->estrip = TArrayI(istrip,svstrip);
1915     };
1916     //
1917     }
1918    
1919     void CaloLevel0::ClearStructs(){
1920     ClearTrkVar();
1921     ClearCommonVar();
1922     }
1923    
1924 mocchiut 1.22 void CaloLevel0::Delete(Option_t *t){
1925     if ( de ) delete de;
1926     delete this;
1927     }
1928    
1929    
1930 mocchiut 1.1 void CaloLevel0::RunClose(){
1931     l0tr->Delete();
1932     ClearStructs();
1933     //
1934     memset(dexy, 0, 2*22*96*sizeof(Float_t));
1935     memset(dexyc, 0, 2*22*96*sizeof(Float_t));
1936     memset(base, 0, 2*22*6*sizeof(Float_t));
1937     memset(sbase, 0, 2*22*6*sizeof(Float_t));
1938 mocchiut 1.9 memset(ctprecor, 0, 2*22*6*sizeof(Float_t));
1939 mocchiut 1.22 memset(ctsicor, 0, 2*22*9*sizeof(Float_t));
1940 mocchiut 1.9 memset(ctneigcor, 0, 2*22*6*sizeof(Float_t));
1941 mocchiut 1.1 //
1942     }
1943    
1944     //
1945     // Private methods
1946     //
1947    
1948     void CaloLevel0::ClearTrkVar(){
1949     clevel2->ncore = 0;
1950     clevel2->qcore = 0.;
1951     clevel2->noint = 0.;
1952     clevel2->ncyl = 0.;
1953     clevel2->qcyl = 0.;
1954     clevel2->qtrack = 0.;
1955     clevel2->qtrackx = 0.;
1956     clevel2->qtracky = 0.;
1957     clevel2->dxtrack = 0.;
1958     clevel2->dytrack = 0.;
1959     clevel2->qlast = 0.;
1960     clevel2->nlast = 0.;
1961     clevel2->qpre = 0.;
1962     clevel2->npre = 0.;
1963     clevel2->qpresh = 0.;
1964     clevel2->npresh = 0.;
1965     clevel2->qlow = 0.;
1966     clevel2->nlow = 0.;
1967     clevel2->qtr = 0.;
1968     clevel2->ntr = 0.;
1969     clevel2->planetot = 0.;
1970     clevel2->qmean = 0.;
1971     clevel2->dX0l = 0.;
1972     clevel2->elen = 0.;
1973     clevel2->selen = 0.;
1974     memset(clevel1->al_p, 0, 5*2*sizeof(Double_t));
1975     memset(clevel2->tibar, 0, 2*22*sizeof(Int_t));
1976     memset(clevel2->tbar, 0, 2*22*sizeof(Float_t));
1977     }
1978    
1979     void CaloLevel0::ClearCommonVar(){
1980     istrip = 0;
1981     clevel2->trigty = -1.;
1982     clevel2->wartrig = 0.;
1983     clevel2->good = 0;
1984     clevel2->nstrip = 0.;
1985 mocchiut 1.13 clevel2->nsatstrip = 0.;
1986 mocchiut 1.1 clevel2->qtot = 0.;
1987     // clevel2->impx = 0.;
1988     // clevel2->impy = 0.;
1989     clevel2->tanx = 0.; // this is correct since it refers to the fortran structure
1990     clevel2->tany = 0.; // this is correct since it refers to the fortran structure
1991     clevel2->qmax = 0.;
1992     clevel2->nx22 = 0.;
1993     clevel2->qx22 = 0.;
1994     memset(clevel2->perr, 0, 4*sizeof(Int_t));
1995     memset(clevel2->swerr, 0, 4*sizeof(Int_t));
1996     memset(clevel2->crc, 0, 4*sizeof(Int_t));
1997     memset(clevel2->qq, 0, 4*sizeof(Int_t));
1998     memset(clevel2->varcfit, 0, 4*sizeof(Float_t));
1999     memset(clevel2->npcfit, 0, 4*sizeof(Int_t));
2000     memset(clevel2->planemax, 0, 2*sizeof(Int_t));
2001 mocchiut 1.6 memset(clevel2->selfdelay, 0, 4*7*sizeof(Int_t));
2002 mocchiut 1.1 memset(clevel2->fmode, 0, 2*sizeof(Int_t));
2003     memset(clevel2->cibar, 0, 2*22*sizeof(Int_t));
2004     memset(clevel2->cbar, 0, 2*22*sizeof(Float_t));
2005     }
2006    
2007     void CaloLevel0::ClearCalibVals(Int_t s){
2008     //
2009     for ( Int_t d=0 ; d<11 ;d++ ){
2010     Int_t pre = -1;
2011     for ( Int_t j=0; j<96 ;j++){
2012     if ( j%16 == 0 ) pre++;
2013     if ( s == 2 ){
2014     calped[0][2*d+1][j] = 0.;
2015     cstwerr[3] = 0.;
2016     cperror[3] = 0.;
2017     calgood[0][2*d+1][j] = 0.;
2018     calthr[0][2*d+1][pre] = 0.;
2019     calrms[0][2*d+1][j] = 0.;
2020     calbase[0][2*d+1][pre] = 0.;
2021     calvar[0][2*d+1][pre] = 0.;
2022     };
2023     if ( s == 3 ){
2024     calped[0][2*d][j] = 0.;
2025     cstwerr[1] = 0.;
2026     cperror[1] = 0.;
2027     calgood[0][2*d][j] = 0.;
2028     calthr[0][2*d][pre] = 0.;
2029     calrms[0][2*d][j] = 0.;
2030     calbase[0][2*d][pre] = 0.;
2031     calvar[0][2*d][pre] = 0.;
2032     };
2033     if ( s == 0 ){
2034     calped[1][2*d][j] = 0.;
2035     cstwerr[0] = 0.;
2036     cperror[0] = 0.;
2037     calgood[1][2*d][j] = 0.;
2038     calthr[1][2*d][pre] = 0.;
2039     calrms[1][2*d][j] = 0.;
2040     calbase[1][2*d][pre] = 0.;
2041     calvar[1][2*d][pre] = 0.;
2042     };
2043     if ( s == 1 ){
2044     calped[1][2*d+1][j] = 0.;
2045     cstwerr[2] = 0.;
2046     cperror[2] = 0.;
2047     calgood[1][2*d+1][j] = 0.;
2048     calthr[1][2*d+1][pre] = 0.;
2049     calrms[1][2*d+1][j] = 0.;
2050     calbase[1][2*d+1][pre] = 0.;
2051     calvar[1][2*d+1][pre] = 0.;
2052     };
2053     };
2054     };
2055     return;
2056     }
2057    
2058 mocchiut 1.7 Int_t CaloLevel0::Update(GL_TABLES *glt, UInt_t atime, Int_t s){
2059 mocchiut 1.1 //
2060 mocchiut 1.7 const TString host = glt->CGetHost();
2061     const TString user = glt->CGetUser();
2062     const TString psw = glt->CGetPsw();
2063     TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
2064     if ( !dbc->IsConnected() ) throw -116;
2065 mocchiut 1.8 stringstream myquery;
2066     myquery.str("");
2067     myquery << "SET time_zone='+0:00'";
2068     dbc->Query(myquery.str().c_str());
2069 mocchiut 1.1 Int_t sgnl = 0;
2070     //
2071     GL_CALO_CALIB *glcalo = new GL_CALO_CALIB();
2072     //
2073     sgnl = 0;
2074     //
2075     idcalib[s] = 0;
2076     fromtime[s] = 0;
2077     totime[s] = 0;
2078     calibno[s] = 0;
2079     ClearCalibVals(s);
2080     //
2081     UInt_t uptime = 0;
2082     //
2083     sgnl = glcalo->Query_GL_CALO_CALIB(atime,uptime,s,dbc);
2084     if ( sgnl < 0 ){
2085     if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
2086     return(sgnl);
2087     };
2088     //
2089     idcalib[s] = glcalo->ID_ROOT_L0;
2090     fromtime[s] = glcalo->FROM_TIME;
2091     if ( glcalo->TO_TIME < atime ){ // calibration is corrupted and we are using the one that preceed the good one
2092     totime[s] = uptime;
2093     } else {
2094     totime[s] = glcalo->TO_TIME;
2095     };
2096     // totime[s] = glcalo->TO_TIME;
2097     calibno[s] = glcalo->EV_ROOT;
2098     //
2099     if ( totime[s] == 0 ){
2100     if ( verbose ) printf(" CALORIMETER - WARNING: data with no associated calibration\n");
2101     ClearCalibVals(s);
2102     sgnl = 100;
2103     };
2104     //
2105     // determine path and name and entry of the calibration file
2106     //
2107     GL_ROOT *glroot = new GL_ROOT();
2108     if ( verbose ) printf("\n");
2109     if ( verbose ) printf(" ** SECTION %i **\n",s);
2110     //
2111     sgnl = glroot->Query_GL_ROOT(idcalib[s],dbc);
2112     if ( sgnl < 0 ){
2113     if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
2114     return(sgnl);
2115     };
2116     //
2117     stringstream name;
2118     name.str("");
2119     name << glroot->PATH.Data() << "/";
2120     name << glroot->NAME.Data();
2121     //
2122     fcalname[s] = (TString)name.str().c_str();
2123     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]);
2124     //
2125     sgnl = LoadCalib(s);
2126     //
2127     if ( sgnl != 0 ) return(sgnl);
2128     delete glcalo;
2129     delete glroot;
2130     //
2131     return(0);
2132     //
2133     }
2134    
2135     Int_t CaloLevel0::LoadCalib(Int_t s){
2136     //
2137     ifstream myfile;
2138     myfile.open(fcalname[s].Data());
2139     if ( !myfile ){
2140     return(-107);
2141     };
2142     myfile.close();
2143     //
2144     TFile *File = new TFile(fcalname[s].Data());
2145     if ( !File ) return(-108);
2146     TTree *tr = (TTree*)File->Get("CalibCalPed");
2147     if ( !tr ) return(-109);
2148     //
2149     TBranch *calo = tr->GetBranch("CalibCalPed");
2150     //
2151     pamela::CalibCalPedEvent *ce = 0;
2152     tr->SetBranchAddress("CalibCalPed", &ce);
2153     //
2154     Long64_t ncalibs = calo->GetEntries();
2155     //
2156     if ( !ncalibs ) return(-110);
2157     //
2158     calo->GetEntry(calibno[s]);
2159     //
2160     if (ce->cstwerr[s] != 0 && ce->cperror[s] == 0 ) {
2161     for ( Int_t d=0 ; d<11 ;d++ ){
2162     Int_t pre = -1;
2163     for ( Int_t j=0; j<96 ;j++){
2164     if ( j%16 == 0 ) pre++;
2165     if ( s == 2 ){
2166     calped[0][2*d+1][j] = ce->calped[3][d][j];
2167     cstwerr[3] = ce->cstwerr[3];
2168     cperror[3] = ce->cperror[3];
2169     calgood[0][2*d+1][j] = ce->calgood[3][d][j];
2170     calthr[0][2*d+1][pre] = ce->calthr[3][d][pre];
2171     calrms[0][2*d+1][j] = ce->calrms[3][d][j];
2172     calbase[0][2*d+1][pre] = ce->calbase[3][d][pre];
2173     calvar[0][2*d+1][pre] = ce->calvar[3][d][pre];
2174     };
2175     if ( s == 3 ){
2176     calped[0][2*d][j] = ce->calped[1][d][j];
2177     cstwerr[1] = ce->cstwerr[1];
2178     cperror[1] = ce->cperror[1];
2179     calgood[0][2*d][j] = ce->calgood[1][d][j];
2180     calthr[0][2*d][pre] = ce->calthr[1][d][pre];
2181     calrms[0][2*d][j] = ce->calrms[1][d][j];
2182     calbase[0][2*d][pre] = ce->calbase[1][d][pre];
2183     calvar[0][2*d][pre] = ce->calvar[1][d][pre];
2184     };
2185     if ( s == 0 ){
2186     calped[1][2*d][j] = ce->calped[0][d][j];
2187     cstwerr[0] = ce->cstwerr[0];
2188     cperror[0] = ce->cperror[0];
2189     calgood[1][2*d][j] = ce->calgood[0][d][j];
2190     calthr[1][2*d][pre] = ce->calthr[0][d][pre];
2191     calrms[1][2*d][j] = ce->calrms[0][d][j];
2192     calbase[1][2*d][pre] = ce->calbase[0][d][pre];
2193     calvar[1][2*d][pre] = ce->calvar[0][d][pre];
2194     };
2195     if ( s == 1 ){
2196     calped[1][2*d+1][j] = ce->calped[2][d][j];
2197     cstwerr[2] = ce->cstwerr[2];
2198     cperror[2] = ce->cperror[2];
2199     calgood[1][2*d+1][j] = ce->calgood[2][d][j];
2200     calthr[1][2*d+1][pre] = ce->calthr[2][d][pre];
2201     calrms[1][2*d+1][j] = ce->calrms[2][d][j];
2202     calbase[1][2*d+1][pre] = ce->calbase[2][d][pre];
2203     calvar[1][2*d+1][pre] = ce->calvar[2][d][pre];
2204     };
2205     };
2206     };
2207     } else {
2208     if ( verbose ) printf(" CALORIMETER - ERROR: problems finding a good calibration in this file! \n\n ");
2209     return(-111);
2210     };
2211     File->Close();
2212     return(0);
2213     }

  ViewVC Help
Powered by ViewVC 1.1.23