/[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.7 - (show annotations) (download)
Sun Sep 9 18:57:26 2007 UTC (17 years, 2 months ago) by mocchiut
Branch: MAIN
CVS Tags: v4r00
Changes since 1.6: +26 -6 lines
Changed core routines to open themselves the DB connection

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

  ViewVC Help
Powered by ViewVC 1.1.23