/[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.27 - (hide annotations) (download)
Wed Nov 25 09:38:30 2009 UTC (15 years ago) by mocchiut
Branch: MAIN
Changes since 1.26: +8 -2 lines
Bug in CaloLevel2::good variable fixed

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

  ViewVC Help
Powered by ViewVC 1.1.23