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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.9 - (show annotations) (download)
Thu Oct 25 13:24:19 2007 UTC (17 years, 4 months ago) by mocchiut
Branch: MAIN
Changes since 1.8: +393 -8 lines
Calorimeter crosstalk code changed in order to use flight calibrations

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

  ViewVC Help
Powered by ViewVC 1.1.23