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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.2 - (show annotations) (download)
Fri Jun 30 09:21:54 2006 UTC (18 years, 5 months ago) by mocchiut
Branch: MAIN
CVS Tags: v1r01beta, v1r01
Changes since 1.1: +10 -0 lines
Memory leak bug fixed

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

  ViewVC Help
Powered by ViewVC 1.1.23