/[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.1.1.1 - (show annotations) (download) (vendor branch)
Fri May 19 13:15:49 2006 UTC (18 years, 6 months ago) by mocchiut
Branch: DarthVader
CVS Tags: v0r01, v0r02, v1r00, start
Changes since 1.1: +0 -0 lines
Imported sources

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

  ViewVC Help
Powered by ViewVC 1.1.23