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