/[PAMELA software]/PamelaLevel2/src/PamLevel2.cpp
ViewVC logotype

Contents of /PamelaLevel2/src/PamLevel2.cpp

Parent Directory Parent Directory | Revision Log Revision Log

Revision 1.105 - (show annotations) (download)
Fri Oct 17 16:22:26 2014 UTC (10 years, 4 months ago) by mocchiut
Branch: MAIN
Changes since 1.104: +1 -1 lines 
ProcessingInfo readout added (not fully implemented yet)
1 #include <PamLevel2.h>
2
3 /////////////////////////////////////////////////////////////////////////////
4 /////////////////////////////////////////////////////////////////////////////
5 /////////////////////////////////////////////////////////////////////////////
6 /////////////////////////////////////////////////////////////////////////////
7
8 void GPamela::Clear() {
9 Irun = 0;
10 Ievnt = 0;
11 Ipa = 0;
12 X0 = 0.;
13 Y0 = 0.;
14 Z0 = 0.;
15 Theta = 0.;
16 Phi = 0.;
17 P0 = 0.;
18 Nthtof = 0;
19 Nthcat = 0;
20 Nthcas = 0;
21 Nthspe = 0;
22 Nstrpx = 0;
23 Nstrpy = 0;
24 Nthcali = 0;
25 Nthcal = 0;
26 Nthnd = 0;
27 Nthcard = 0;
28 }
29
30 void GPamela::Delete() {
31
32 Clear();
33 /* EM all these are in the stack not in the heap, so no need to delete by hand...
34 delete[] Ipltof;
35 delete[] Ipaddle;
36 delete[] Ipartof;
37 delete[] Xintof;
38 delete[] Yintof;
39 delete[] Zintof;
40 delete[] Xouttof;
41 delete[] Youttof;
42 delete[] Zouttof;
43 delete[] Ereltof;
44 delete[] Timetof;
45 delete[] Pathtof;
46 delete[] P0tof;
47 delete[] Iparcat;
48 delete[] Icat;
49 delete[] Xincat;
50 delete[] Yincat;
51 delete[] Zincat;
52 delete[] Xoutcat;
53 delete[] Youtcat;
54 delete[] Zoutcat;
55 delete[] Erelcat;
56 delete[] Timecat;
57 delete[] Pathcat;
58 delete[] P0cat;
59 delete[] Iparcas;
60 delete[] Icas;
61 delete[] Xincas;
62 delete[] Yincas;
63 delete[] Zincas;
64 delete[] Xoutcas;
65 delete[] Youtcas;
66 delete[] Zoutcas;
67 delete[] Erelcas;
68 delete[] Timecas;
69 delete[] Pathcas;
70 delete[] P0cas;
71 delete[] Iparspe;
72 delete[] Itrpb;
73 delete[] Itrsl;
74 delete[] Itspa;
75 delete[] Xinspe;
76 delete[] Yinspe;
77 delete[] Zinspe;
78 delete[] Xoutspe;
79 delete[] Youtspe;
80 delete[] Zoutspe;
81 delete[] Xavspe;
82 delete[] Yavspe;
83 delete[] Zavspe;
84 delete[] Erelspe;
85 delete[] Pathspe;
86 delete[] P0spe;
87 delete[] Nxmult;
88 delete[] Nymult;
89 delete[] Npstripx;
90 delete[] Ntstripx;
91 delete[] Istripx;
92 delete[] Qstripx;
93 delete[] Xstripx;
94 delete[] Npstripy;
95 delete[] Ntstripy;
96 delete[] Istripy;
97 delete[] Qstripy;
98 delete[] Ystripy;
99 delete[] Icaplane;
100 delete[] Icastrip;
101 delete[] Icamod;
102 delete[] Enestrip;
103 delete[] Icapl;
104 delete[] Icasi;
105 delete[] Icast;
106 delete[] Xincal;
107 delete[] Yincal;
108 delete[] Zincal;
109 delete[] Erelcal;
110 delete[] Itubend;
111 delete[] Iparnd;
112 delete[] Xinnd;
113 delete[] Yinnd;
114 delete[] Zinnd;
115 delete[] Xoutnd;
116 delete[] Youtnd;
117 delete[] Zoutnd;
118 delete[] Erelnd;
119 delete[] Timend;
120 delete[] Pathnd;
121 delete[] P0nd;
122 delete[] Iparcard;
123 delete[] Icard;
124 delete[] Xincard;
125 delete[] Yincard;
126 delete[] Zincard;
127 delete[] Xoutcard;
128 delete[] Youtcard;
129 delete[] Zoutcard;
130 delete[] Erelcard;
131 delete[] Timecard;
132 delete[] Pathcard;
133 delete[] P0card;
134 */
135 }
136
137
138 void GPamela::SetBranchAddress(TChain* fhBookTree) {
139
140 // cout << "fhBookTree "<<fhBookTree << endl;
141 // prepare tree
142 fhBookTree->SetBranchAddress("Irun", &Irun);
143 fhBookTree->SetBranchAddress("Ievnt", &Ievnt);
144 fhBookTree->SetBranchAddress("Ipa", &Ipa);
145 fhBookTree->SetBranchAddress("X0", &X0);
146 fhBookTree->SetBranchAddress("Y0", &Y0);
147 fhBookTree->SetBranchAddress("Z0", &Z0);
148 fhBookTree->SetBranchAddress("Theta", &Theta);
149 fhBookTree->SetBranchAddress("Phi", &Phi);
150 fhBookTree->SetBranchAddress("P0", &P0);
151 fhBookTree->SetBranchAddress("Nthtof", &Nthtof);
152 fhBookTree->SetBranchAddress("Ipltof", Ipltof);
153 fhBookTree->SetBranchAddress("Ipaddle", Ipaddle);
154 fhBookTree->SetBranchAddress("Ipartof", Ipartof);
155 fhBookTree->SetBranchAddress("Xintof", Xintof);
156 fhBookTree->SetBranchAddress("Yintof", Yintof);
157 fhBookTree->SetBranchAddress("Zintof", Zintof);
158 fhBookTree->SetBranchAddress("Xouttof", Xouttof);
159 fhBookTree->SetBranchAddress("Youttof", Youttof);
160 fhBookTree->SetBranchAddress("Zouttof", Zouttof);
161 fhBookTree->SetBranchAddress("Ereltof", Ereltof);
162 fhBookTree->SetBranchAddress("Timetof", Timetof);
163 fhBookTree->SetBranchAddress("Pathtof", Pathtof);
164 fhBookTree->SetBranchAddress("P0tof", P0tof);
165 fhBookTree->SetBranchAddress("Nthcat", &Nthcat);
166 fhBookTree->SetBranchAddress("Iparcat", Iparcat);
167 fhBookTree->SetBranchAddress("Icat", Icat);
168 fhBookTree->SetBranchAddress("Xincat", Xincat);
169 fhBookTree->SetBranchAddress("Yincat", Yincat);
170 fhBookTree->SetBranchAddress("Zincat", Zincat);
171 fhBookTree->SetBranchAddress("Xoutcat", Xoutcat);
172 fhBookTree->SetBranchAddress("Youtcat", Youtcat);
173 fhBookTree->SetBranchAddress("Zoutcat", Zoutcat);
174 fhBookTree->SetBranchAddress("Erelcat", Erelcat);
175 fhBookTree->SetBranchAddress("Timecat", Timecat);
176 fhBookTree->SetBranchAddress("Pathcat", Pathcat);
177 fhBookTree->SetBranchAddress("P0cat", P0cat);
178 fhBookTree->SetBranchAddress("Nthcas", &Nthcas);
179 fhBookTree->SetBranchAddress("Iparcas", Iparcas);
180 fhBookTree->SetBranchAddress("Icas", Icas);
181 fhBookTree->SetBranchAddress("Xincas", Xincas);
182 fhBookTree->SetBranchAddress("Yincas", Yincas);
183 fhBookTree->SetBranchAddress("Zincas", Zincas);
184 fhBookTree->SetBranchAddress("Xoutcas", Xoutcas);
185 fhBookTree->SetBranchAddress("Youtcas", Youtcas);
186 fhBookTree->SetBranchAddress("Zoutcas", Zoutcas);
187 fhBookTree->SetBranchAddress("Erelcas", Erelcas);
188 fhBookTree->SetBranchAddress("Timecas", Timecas);
189 fhBookTree->SetBranchAddress("Pathcas", Pathcas);
190 fhBookTree->SetBranchAddress("P0cas", P0cas);
191 fhBookTree->SetBranchAddress("Nthspe", &Nthspe);
192 fhBookTree->SetBranchAddress("Iparspe", Iparspe);
193 fhBookTree->SetBranchAddress("Itrpb", Itrpb);
194 fhBookTree->SetBranchAddress("Itrsl", Itrsl);
195 fhBookTree->SetBranchAddress("Itspa", Itspa);
196 fhBookTree->SetBranchAddress("Xinspe", Xinspe);
197 fhBookTree->SetBranchAddress("Yinspe", Yinspe);
198 fhBookTree->SetBranchAddress("Zinspe", Zinspe);
199 fhBookTree->SetBranchAddress("Xoutspe", Xoutspe);
200 fhBookTree->SetBranchAddress("Youtspe", Youtspe);
201 fhBookTree->SetBranchAddress("Zoutspe", Zoutspe);
202 fhBookTree->SetBranchAddress("Xavspe", Xavspe);
203 fhBookTree->SetBranchAddress("Yavspe", Yavspe);
204 fhBookTree->SetBranchAddress("Zavspe", Zavspe);
205 fhBookTree->SetBranchAddress("Erelspe", Erelspe);
206 fhBookTree->SetBranchAddress("Pathspe", Pathspe);
207 fhBookTree->SetBranchAddress("P0spe", P0spe);
208 fhBookTree->SetBranchAddress("Nxmult", Nxmult);
209 fhBookTree->SetBranchAddress("Nymult", Nymult);
210 fhBookTree->SetBranchAddress("Nstrpx", &Nstrpx);
211 fhBookTree->SetBranchAddress("Npstripx", Npstripx);
212 fhBookTree->SetBranchAddress("Ntstripx", Ntstripx);
213 fhBookTree->SetBranchAddress("Istripx", Istripx);
214 fhBookTree->SetBranchAddress("Qstripx", Qstripx);
215 fhBookTree->SetBranchAddress("Xstripx", Xstripx);
216 fhBookTree->SetBranchAddress("Nstrpy", &Nstrpy);
217 fhBookTree->SetBranchAddress("Npstripy", Npstripy);
218 fhBookTree->SetBranchAddress("Ntstripy", Ntstripy);
219 fhBookTree->SetBranchAddress("Istripy", Istripy);
220 fhBookTree->SetBranchAddress("Qstripy", Qstripy);
221 fhBookTree->SetBranchAddress("Ystripy", Ystripy);
222 fhBookTree->SetBranchAddress("Nthcali", &Nthcali);
223 fhBookTree->SetBranchAddress("Icaplane", Icaplane);
224 fhBookTree->SetBranchAddress("Icastrip", Icastrip);
225 fhBookTree->SetBranchAddress("Icamod", Icamod);
226 fhBookTree->SetBranchAddress("Enestrip", Enestrip);
227 fhBookTree->SetBranchAddress("Nthcal", &Nthcal);
228 fhBookTree->SetBranchAddress("Icapl", Icapl);
229 fhBookTree->SetBranchAddress("Icasi", Icasi);
230 fhBookTree->SetBranchAddress("Icast", Icast);
231 fhBookTree->SetBranchAddress("Xincal", Xincal);
232 fhBookTree->SetBranchAddress("Yincal", Yincal);
233 fhBookTree->SetBranchAddress("Zincal", Zincal);
234 fhBookTree->SetBranchAddress("Erelcal", Erelcal);
235 fhBookTree->SetBranchAddress("Nthnd", &Nthnd);
236 fhBookTree->SetBranchAddress("Itubend", Itubend);
237 fhBookTree->SetBranchAddress("Iparnd", Iparnd);
238 fhBookTree->SetBranchAddress("Xinnd", Xinnd);
239 fhBookTree->SetBranchAddress("Yinnd", Yinnd);
240 fhBookTree->SetBranchAddress("Zinnd", Zinnd);
241 fhBookTree->SetBranchAddress("Xoutnd", Xoutnd);
242 fhBookTree->SetBranchAddress("Youtnd", Youtnd);
243 fhBookTree->SetBranchAddress("Zoutnd", Zoutnd);
244 fhBookTree->SetBranchAddress("Erelnd", Erelnd);
245 fhBookTree->SetBranchAddress("Timend", Timend);
246 fhBookTree->SetBranchAddress("Pathnd", Pathnd);
247 fhBookTree->SetBranchAddress("P0nd", P0nd);
248 fhBookTree->SetBranchAddress("Nthcard", &Nthcard);
249 fhBookTree->SetBranchAddress("Iparcard", Iparcard);
250 fhBookTree->SetBranchAddress("Icard", Icard);
251 fhBookTree->SetBranchAddress("Xincard", Xincard);
252 fhBookTree->SetBranchAddress("Yincard", Yincard);
253 fhBookTree->SetBranchAddress("Zincard", Zincard);
254 fhBookTree->SetBranchAddress("Xoutcard", Xoutcard);
255 fhBookTree->SetBranchAddress("Youtcard", Youtcard);
256 fhBookTree->SetBranchAddress("Zoutcard", Zoutcard);
257 fhBookTree->SetBranchAddress("Erelcard", Erelcard);
258 fhBookTree->SetBranchAddress("Timecard", Timecard);
259 fhBookTree->SetBranchAddress("Pathcard", Pathcard);
260 fhBookTree->SetBranchAddress("P0card", P0card);
261
262 // fhBookTree->SetBranchStatus("*",0);
263
264 }
265
266 ClassImp( GPamela);
267
268 /////////////////////////////////////////////////////////////////////////////
269 /////////////////////////////////////////////////////////////////////////////
270 /////////////////////////////////////////////////////////////////////////////
271 /////////////////////////////////////////////////////////////////////////////
272 //--------------------------------------
273 //
274 //
275 //--------------------------------------
276 /**
277 * Default constructor
278 */
279 PamTrack::PamTrack() {
280 trk_track = 0;
281 trk_ext_track = 0;
282 calo_track = 0;
283 tof_track = 0;
284 orb_track = 0;
285 candeleteobj = 0;
286 pscore = 0;
287 iscore = 0;
288 }
289 ;
290 //--------------------------------------
291 //
292 //
293 //--------------------------------------
294 /**
295 * Constructor
296 */
297 PamTrack::PamTrack(TrkTrack* t, CaloTrkVar* c, ToFTrkVar* o, OrbitalInfoTrkVar *r) {
298
299 trk_ext_track = 0;
300 trk_track = 0;
301 calo_track = 0;
302 tof_track = 0;
303 orb_track = 0;
304 // if(t)trk_track = new TrkTrack(*t);
305 // if(c)calo_track = new CaloTrkVar(*c);
306 // if(o)tof_track = new ToFTrkVar(*o);
307 // if (t)
308 // trk_track = t;
309 // if (c)
310 // calo_track = c;
311 // if (o)
312 // tof_track = o;
313 // if (r)
314 // orb_track = r;
315
316 // candeleteobj = 0;
317
318
319 if (t){
320 trk_track = new TrkTrack(*t);
321 trk_ext_track = new ExtTrack(*t);//NB!! ha dimensione 6 invece che 8
322 }
323 if (c)
324 calo_track = new CaloTrkVar(*c);
325 if (o)
326 tof_track = new ToFTrkVar(*o);
327 if (r)
328 orb_track = new OrbitalInfoTrkVar(*r);
329 candeleteobj = 1;
330
331 }
332 ;
333 /**
334 * Constructor
335 */
336 PamTrack::PamTrack(ExtTrack* t, CaloTrkVar* c, ToFTrkVar* o, OrbitalInfoTrkVar *r) {
337
338 trk_ext_track = 0;
339 trk_track = 0;
340 calo_track = 0;
341 tof_track = 0;
342 orb_track = 0;
343 // if(t)trk_track = new TrkTrack(*t);
344 // if(c)calo_track = new CaloTrkVar(*c);
345 // if(o)tof_track = new ToFTrkVar(*o);
346 // if (t)
347 // trk_track = t;
348 // if (c)
349 // calo_track = c;
350 // if (o)
351 // tof_track = o;
352 // if (r)
353 // orb_track = r;
354
355 // candeleteobj = 0;
356
357
358 if (t){
359 //// trk_track = new TrkTrack(*t);//in this case TrkTrack object remains null
360 trk_ext_track = new ExtTrack(*t);
361 }
362 if (c)
363 calo_track = new CaloTrkVar(*c);
364 if (o)
365 tof_track = new ToFTrkVar(*o);
366 if (r)
367 orb_track = new OrbitalInfoTrkVar(*r);
368 candeleteobj = 1;
369 pscore = 0;
370 iscore = 0;
371
372 }
373 ;
374 PamTrack::PamTrack(const PamTrack& track) {
375
376 TrkTrack *t = track.trk_track;
377 ExtTrack *e = track.trk_ext_track;
378 CaloTrkVar *c = track.calo_track;
379 ToFTrkVar *o = track.tof_track;
380 OrbitalInfoTrkVar *r = track.orb_track;
381
382 trk_ext_track = 0;
383 trk_track = 0;
384 calo_track = 0;
385 tof_track = 0;
386 orb_track = 0;
387 if(e)
388 trk_ext_track = new ExtTrack(*e);
389 if (t)
390 trk_track = new TrkTrack(*t);
391 if (c)
392 calo_track = new CaloTrkVar(*c);
393 if (o)
394 tof_track = new ToFTrkVar(*o);
395 if (r)
396 orb_track = new OrbitalInfoTrkVar(*r);
397
398 candeleteobj = 1;
399 pscore = 0;
400 iscore = 0;
401
402 }
403 void PamTrack::Clear(Option_t *option) {
404
405 // cout << "PamTrack::Clear( "<<option<<" ) "<<candeleteobj<<endl;
406 if (candeleteobj) {
407
408 if (trk_ext_track)
409 trk_ext_track->ExtTrack::Clear(option);
410 if (trk_track)
411 trk_track->TrkTrack::Clear();
412 if (calo_track)
413 calo_track->CaloTrkVar::Clear();//???
414 if (tof_track)
415 tof_track->ToFTrkVar::Clear();//???
416 if (orb_track)
417 orb_track->OrbitalInfoTrkVar::Clear();//???
418 }
419 else {
420 trk_ext_track = 0;
421 trk_track = 0;
422 calo_track = 0;
423 tof_track = 0;
424 orb_track = 0;
425 }
426 pscore = 0;
427 iscore = 0;
428
429 }
430 void PamTrack::Delete() {
431 // cout << "PamTrack::Delete() "<<candeleteobj<<endl;
432 if (candeleteobj) {
433 if (trk_ext_track) {
434 trk_ext_track->ExtTrack::Clear("C");
435 delete trk_ext_track;
436 }
437 if (trk_track) {
438 trk_track->TrkTrack::Clear();
439 delete trk_track;
440 }
441 if (calo_track) {
442 calo_track->CaloTrkVar::Clear();//???
443 delete calo_track;
444 }
445 if (tof_track) {
446 tof_track->ToFTrkVar::Clear();//???
447 delete tof_track;
448 }
449 if (orb_track) {
450 orb_track->OrbitalInfoTrkVar::Clear();//???
451 delete orb_track;
452 }
453 }
454 else {
455 Clear();
456 }
457 }
458
459
460
461
462 //--------------------------------------
463 //
464 //
465 //--------------------------------------
466 /**
467 * Default Constructor
468 */
469 PamLevel2::PamLevel2() {
470 Initialize();
471 }
472
473
474 /**
475 * Constructor
476 * @param ddir Name of directory where level2 files are stored.
477 * @param list Name of an ascii file containing the list of file names
478 * @param detlist Options to chose what to load.
479 * Possible options are:
480 * +AUTO --> load all trees/branches in the input files
481 * +(-)ALL --> inlcude(exclude) all trees/branches
482 * +(-)TRK1+(-)TRK2+(-)CAL1+(-)CAL2+(-)TOF+(-)TRG+(-)ND+(-)S4+(-)ORB+(-)AC --> inlcude(exclude) trees and branches
483 * +(-)TRK0 --> include(exclude) tracker level0 tree
484 * +(-)GP --> include exclude GPAMELA output tree
485 * If no options are specified, the default is assumed. Default is:
486 * +TRK2+CAL2+CAL1+TOF+TRG+ND+AC+S4+ORB
487 */
488 PamLevel2::PamLevel2(TString ddir, TString llist, TString detlist) {
489 Initialize();
490 TList* listf = GetListOfLevel2Files(ddir, llist);
491 if (listf)
492 GetPamTree(listf, detlist);
493 if (listf)
494 GetRunTree(listf);
495 SetMaxShift(-1);
496 }
497
498
499 PamLevel2::PamLevel2(TString ddir, TList *llist, TString detlist) {
500 Initialize();
501 GetPamTree(llist, detlist);
502 GetRunTree(llist);
503 SetMaxShift(-1);
504 }
505
506 /**
507 * Constructor
508 * @param ddir Name of directory where level2 files are stored.
509 * @param list Name of an ascii file containing the list of file names
510 * Default trees/branches are loaded. Default is:
511 * +TRK2+CAL2+CAL1+TOF+TRG+ND+AC+S4+ORB
512 */
513 PamLevel2::PamLevel2(TString ddir, TString llist) {
514 Initialize();
515 TList* listf = GetListOfLevel2Files(ddir, llist);
516 cout << "GetPamTree: "<<endl;
517 GetPamTree(listf, "");
518 cout << "GetRunTree: "<<endl;
519 GetRunTree(listf);
520 SetMaxShift(-1);
521 }
522
523
524 void PamLevel2::Initialize() {
525
526 h0_obj = 0;
527 trk0_obj = 0;
528 calo0_obj = 0;
529
530 trk2_obj = 0;
531 trk1_obj = 0;
532 trkh_obj = 0;
533 calo1_obj = 0;
534 calo2_obj = 0;
535 tof2_obj = 0;
536 trig_obj = 0;
537 s4_obj = 0;
538 nd_obj = 0;
539 ac_obj = 0;
540 orb2_obj = 0;
541 gp_obj = 0;
542
543 extAlgFlag=0;
544
545 trk_ext_obj = 0;
546 trk_ext_nuc_obj = 0;
547 trk_nuc_obj = 0;
548
549 calo_ext_obj = 0;
550 calo_ext_nuc_obj = 0;
551 calo_nuc_obj = 0;
552
553 tof_ext_obj = 0;
554 tof_ext_nuc_obj = 0;
555 tof_nuc_obj = 0;
556
557 orb_ext_obj = 0;
558 orb_ext_nuc_obj = 0;
559 orb_nuc_obj = 0;
560
561 trk2_nuc_obj = 0;
562 calo2_nuc_obj = 0;
563 tof2_nuc_obj = 0;
564 orb2_nuc_obj = 0;
565
566
567 run_obj = 0;//new GL_RUN();
568 soft_obj = 0;// Emiliano
569 proc_obj = 0;// Emiliano
570 irun = -1LL;
571 irunt = -1LL;
572 totrunentry = 0LL;
573 totrunentrymax = 0LL;
574 totrunentrymin = 0LL;
575 runfirstentry = 0LL;
576 runlastentry = 0LL;
577 gltsync = 0; // Emiliano
578 fUpdateRunInfo = true; // Emiliano
579 fUseDBinRunInfo = true; // Emiliano
580 fDiscarded = false; //EM
581 isSync = false; // by default assume that the level2 file(s) is(are) not sinchronized between L0/DB and L2, that is we miss some packets in L2 due to nested/DV-skipped events
582 il0entry = 0LL;
583 // hasL0EE = true;
584
585 l0_file = NULL;
586 l0_tree = NULL;
587 iroot = -1;
588 dbc = 0;
589
590 prevshift = 0;
591 yprevshift = 0;
592 maxshift = 10; //EMILIANO now overridden by SetMaxShift(-1) called by constructors
593
594 run_tree = NULL;
595 run_tree_clone = NULL;
596
597 proc_tree = NULL;
598 proc_tree_clone = NULL;
599
600 sel_tree = NULL;
601 sel_tree_clone = NULL;
602
603 irunentry = -1LL;
604 pam_tree = NULL;
605 for (Int_t i = 0; i < NCLONES; i++)
606 pam_tree_clone[i] = NULL;
607
608 totdltime[0] = 0LL;
609 totdltime[1] = 0LL;
610 totdltime[2] = 0LL;
611
612 host = "mysql://localhost/pamelaprod";
613 user = "anonymous";
614 psw = "";
615 const char *pamdbhost = gSystem->Getenv("PAM_DBHOST");
616 const char *pamdbuser = gSystem->Getenv("PAM_DBUSER");
617 const char *pamdbpsw = gSystem->Getenv("PAM_DBPSW");
618 if (!pamdbhost)
619 pamdbhost = "";
620 if (!pamdbuser)
621 pamdbuser = "";
622 if (!pamdbpsw)
623 pamdbpsw = "";
624 if (strcmp(pamdbhost, ""))
625 host = pamdbhost;
626 if (strcmp(pamdbuser, ""))
627 user = pamdbuser;
628 if (strcmp(pamdbpsw, ""))
629 psw = pamdbpsw;
630
631 // sorted_tracks = 0;//new TRefArray();
632
633 CAL0 = false;
634 CAL1 = true;
635 CAL2 = true;
636 TRK2 = true;
637 TRK1 = false;
638 TRK0 = false;
639 TRKh = false;
640 TRG = true;
641 TOF = true;
642 TOF0 = false;
643 S4 = true;
644 ND = true;
645 AC = true;
646 ORB = true;
647 PROC = true;
648 GP = false;
649
650 EXT = false;
651 NUC = false;
652 trkAlg = "STD";//default tracking algorythm
653
654 RUN = true;
655
656 SELLI = -1;
657
658 ISGP = false;
659
660 DBG = false;
661
662 tsorted = 0;
663 timage = 0;
664 text = 0 ;
665
666 tsorted_nuc = 0;
667 timage_nuc = 0;
668 text_nuc = 0 ;
669
670 howtosort = "+CAL+TOF";
671 //howtosort = "+TOF";
672 sortthr = 100.;
673
674 issorted = false;
675 lastsorted = -1;
676 issorted_new = false;
677 lastsorted_new = -1;
678
679 }
680 ;
681 /**
682 * Delete the event (act as Dtor)
683 */
684 void PamLevel2::Delete() {
685
686 if (run_obj)
687 delete run_obj;
688 if (soft_obj)
689 delete soft_obj; //Emiliano
690 if (proc_obj)
691 delete proc_obj; //Emiliano
692
693 // cout << "void PamLevel2::Clear()"<<endl;
694 if (h0_obj)
695 delete h0_obj;
696 if (trk0_obj)
697 delete trk0_obj;
698 if (calo0_obj)
699 delete calo0_obj;
700 if (trk1_obj)
701 delete trk1_obj;
702 if (trk2_obj)
703 delete trk2_obj;
704 if (trkh_obj)
705 delete trkh_obj;
706 if (calo1_obj)
707 delete calo1_obj;
708 if (calo2_obj)
709 delete calo2_obj;
710 if (tof2_obj)
711 delete tof2_obj;
712 if (trig_obj)
713 delete trig_obj;
714 if (s4_obj)
715 delete s4_obj;
716 if (nd_obj)
717 delete nd_obj;
718 if (ac_obj)
719 delete ac_obj;
720 if (orb2_obj)
721 delete orb2_obj;
722 if (gp_obj)
723 delete gp_obj;
724
725 if(trk_nuc_obj)trk_nuc_obj->Delete();
726 if(trk_ext_obj)trk_ext_obj->Delete();
727 if(trk_ext_nuc_obj)trk_ext_nuc_obj->Delete();
728
729 if(calo_nuc_obj)calo_nuc_obj->Delete();
730 if(calo_ext_obj)calo_ext_obj->Delete();
731 if(calo_ext_nuc_obj)calo_ext_nuc_obj->Delete();
732
733 if(tof_nuc_obj)tof_nuc_obj->Delete();
734 if(tof_ext_obj)tof_ext_obj->Delete();
735 if(tof_ext_nuc_obj)tof_ext_nuc_obj->Delete();
736
737 if(orb_nuc_obj)orb_nuc_obj->Delete();
738 if(orb_ext_obj)orb_ext_obj->Delete();
739 if(orb_ext_nuc_obj)orb_ext_nuc_obj->Delete();
740
741
742 if(trk2_nuc_obj)trk2_nuc_obj->Delete();;
743 if( calo2_nuc_obj)calo2_nuc_obj->Delete();;
744 if(tof2_nuc_obj)tof2_nuc_obj->Delete();;
745 if(orb2_nuc_obj)orb2_nuc_obj->Delete();;
746
747
748
749 if (tsorted) {
750 tsorted->Delete();
751 delete tsorted;
752 }
753 if (timage) {
754 timage->Delete();
755 delete timage;
756 }
757 if (text) {
758 text->Delete();
759 delete text;
760 }
761 if (tsorted_nuc) {
762 tsorted_nuc->Delete();
763 delete tsorted_nuc;
764 }
765 if (timage_nuc) {
766 timage_nuc->Delete();
767 delete timage_nuc;
768 }
769 if (text_nuc) {
770 text_nuc->Delete();
771 delete text_nuc;
772 }
773
774
775
776 if (dbc) {
777 dbc->Close();
778 delete dbc;
779 dbc=0;
780 }
781
782 if (gltsync)
783 delete gltsync;
784
785 if (l0_file)
786 l0_file->Close();
787 // if(pam_tree)pam_tree->Delete();;
788
789 if (pam_tree) {
790 //
791 // we have first to find which chains we have to delete, then delete the main chain and only after delete the friend chains. Any other order causes a segfault...
792 //
793 TList *temp = pam_tree->GetListOfFriends();
794 TList *contents = new TList; // create chain friend list
795 contents->SetOwner();
796 TIter next(temp);
797 TChain *questo = 0;
798 while ((questo = (TChain*) next())) {
799 TString name = questo->GetName();
800 contents->Add((TChain*) gROOT->FindObject(name.Data()));// add object to the list
801 };
802 //
803 // deleting the main chain
804 //
805 pam_tree->Delete();
806 //
807 // deleting the friends...
808 //
809 TIter next2(contents);
810 TChain *questa = 0;
811 while ( (questa = (TChain*)next2()) ) {
812 TString name = questa->GetName();
813 questa->Delete();
814 questa = NULL;
815 };
816 //
817 };
818 pam_tree = NULL;
819
820 if (run_tree)
821 run_tree->Delete();;
822 if (sel_tree)
823 sel_tree->Delete();;
824
825 // The following lines are commented out since they may generate a double delete error
826 // if the file containing the clone trees is closed. This is because the file owns the
827 // clone trees which are written into it, so it will delete them when it is closed; if
828 // also PamLevel2 will try to delete these trees, a double delete error will be generated
829 // when exiting from analysis program. (Nicola 28/11/2011)
830
831 /*for (Int_t i = 0; i < NCLONES; i++)
832 if (pam_tree_clone[i])
833 pam_tree_clone[i]->Delete();;
834 if (run_tree_clone)
835 run_tree_clone->Delete();;
836 if (sel_tree_clone)
837 sel_tree_clone->Delete();;*/
838
839 if (irunoffset)
840 delete[] irunoffset;
841
842
843 Initialize();
844
845 }
846 ;
847
848 /**
849 * Clear the event (NB! does not deallocate objects)
850 */
851 void PamLevel2::Clear() {
852
853 // cout << "void PamLevel2::Clear()"<<endl;
854
855 //
856 // This method is called once for every entry but RunInfo and SoftInfo do not change until the next run so we cannot clear them here unless we don't
857 // want to load them for each event even if they are the same...
858 //
859 // if(run_obj)delete run_obj;
860 // if(run_obj) run_obj->Clear(); // Emiliano: Do not deallocate run_obj here, it will give segmentation fault! call clear instead
861 // if(soft_obj) soft_obj->Clear();
862
863 if (h0_obj)
864 h0_obj->Clear();
865 // if(trk0_obj) trk0_obj->Clear();
866 if (trk1_obj)
867 trk1_obj->Clear();
868 if (trk2_obj)
869 trk2_obj->Clear();
870 if (trkh_obj)
871 trkh_obj->Clear();
872 if (calo0_obj)
873 calo0_obj->Clear();
874 if (calo1_obj)
875 calo1_obj->Clear();
876 if (calo2_obj)
877 calo2_obj->Clear();
878 if (tof2_obj)
879 tof2_obj->Clear();
880 if (trig_obj)
881 trig_obj->Clear();
882 if (s4_obj)
883 s4_obj->Clear();
884 if (nd_obj)
885 nd_obj->Clear();
886 if (ac_obj)
887 ac_obj->Clear();
888 if (orb2_obj)
889 orb2_obj->Clear();
890 if (gp_obj)
891 gp_obj->Clear();
892 if (proc_obj)
893 proc_obj->Clear();
894
895 // if(sorted_tracks)sorted_tracks->Clear();
896 // sorted_tracks.Clear();
897
898 if(trk_nuc_obj)trk_nuc_obj->Clear();
899 if(trk_ext_obj)trk_ext_obj->Clear();
900 if(trk_ext_nuc_obj)trk_ext_nuc_obj->Clear();
901
902 if(calo_nuc_obj)calo_nuc_obj->Clear();
903 if(calo_ext_obj)calo_ext_obj->Clear();
904 if(calo_ext_nuc_obj)calo_ext_nuc_obj->Clear();
905
906 if(tof_nuc_obj)tof_nuc_obj->Clear();
907 if(tof_ext_obj)tof_ext_obj->Clear();
908 if(tof_ext_nuc_obj)tof_ext_nuc_obj->Clear();
909
910 if(orb_nuc_obj)orb_nuc_obj->Clear();
911 if(orb_ext_obj)orb_ext_obj->Clear();
912 if(orb_ext_nuc_obj)orb_ext_nuc_obj->Clear();
913
914 if(trk2_nuc_obj)trk2_nuc_obj->Clear();;
915 if( calo2_nuc_obj)calo2_nuc_obj->Clear();;
916 if(tof2_nuc_obj)tof2_nuc_obj->Clear();;
917 if(orb2_nuc_obj)orb2_nuc_obj->Clear();;
918
919 if (tsorted)tsorted->Delete();
920 if (timage)timage->Delete();
921 if (text) text->Delete();
922
923 if (tsorted_nuc)tsorted_nuc->Delete();
924 if (timage_nuc)timage_nuc->Delete();
925 if (text_nuc) text_nuc->Delete();
926 }
927 ;
928
929 void PamLevel2::Reset() {
930 //
931 // First of all clear everything
932 //
933 Clear();
934 //
935 // close and reset chains and pointers
936 //
937 if (pam_tree) {
938 //
939 // we have first to find which chains we have to delete, then delete the main chain and only after delete the friend chains. Any other order causes a segfault...
940 //
941 TList *temp = pam_tree->GetListOfFriends();
942 TList *contents = new TList; // create chain friend list
943 contents->SetOwner();
944 TIter next(temp);
945 TChain *questo = 0;
946 while ((questo = (TChain*) next())) {
947 TString name = questo->GetName();
948 contents->Add((TChain*) gROOT->FindObject(name.Data()));// add object to the list
949 };
950 //
951 // deleting the main chain
952 //
953 pam_tree->Delete();
954 //
955 // deleting the friends...
956 //
957 TIter next2(contents);
958 TChain *questa = 0;
959 while ( (questa = (TChain*) next2()) ) {
960 TString name = questa->GetName();
961 questa->Delete();
962 questa = NULL;
963 };
964 //
965 };
966 pam_tree = NULL;
967 //
968 if (run_tree)
969 run_tree->Delete();;
970 run_tree = NULL;
971 if (sel_tree)
972 sel_tree->Delete();;
973 sel_tree = NULL;
974
975 if (proc_tree)
976 proc_tree->Delete();
977 proc_tree = NULL;
978 //
979 // Close file
980 //
981 if (l0_file)
982 l0_file->Close("R");
983 l0_file = NULL;
984 //
985 h0_obj = 0;
986 trk0_obj = 0;
987 calo0_obj = 0;
988 //
989 trk2_obj = 0;
990 trk1_obj = 0;
991 trkh_obj = 0;
992 calo1_obj = 0;
993 calo2_obj = 0;
994 tof2_obj = 0;
995 trig_obj = 0;
996 s4_obj = 0;
997 nd_obj = 0;
998 ac_obj = 0;
999 orb2_obj = 0;
1000 gp_obj = 0;
1001 proc_obj = 0;
1002
1003 trk_ext_obj = 0;
1004 trk_ext_nuc_obj = 0;
1005 trk_nuc_obj = 0;
1006
1007 calo_ext_obj = 0;
1008 calo_ext_nuc_obj = 0;
1009 calo_nuc_obj = 0;
1010
1011 tof_ext_obj = 0;
1012 tof_ext_nuc_obj = 0;
1013 tof_nuc_obj = 0;
1014
1015 orb_ext_obj = 0;
1016 orb_ext_nuc_obj = 0;
1017 orb_nuc_obj = 0;
1018
1019 trk2_nuc_obj = 0;
1020 calo2_nuc_obj = 0;
1021 tof2_nuc_obj = 0;
1022 orb2_nuc_obj = 0;
1023
1024 trk2_nuc_obj = 0;
1025 calo2_nuc_obj = 0;
1026 tof2_nuc_obj = 0;
1027 orb2_nuc_obj = 0;
1028 //
1029 // Reset run pointers
1030 //
1031 run_obj = 0;//new GL_RUN();
1032 soft_obj = 0;// Emiliano
1033 proc_obj = 0;// Emiliano
1034 irun = -1;
1035 irunt = -1;
1036 totrunentry = 0LL;
1037 totrunentrymax = 0LL;
1038 totrunentrymin = 0LL;
1039 runfirstentry = 0ULL;
1040 runlastentry = 0ULL;
1041 prevabstime = 0ULL;
1042 prevpktnum = 0;
1043 abstime = 0ULL;
1044 pktnum = 0;
1045 isFragment = false;
1046 //
1047 totdltime[0] = 0LL;
1048 totdltime[1] = 0LL;
1049 totdltime[2] = 0LL;
1050 //
1051 }
1052 ;
1053
1054 Bool_t PamLevel2::IsGood(Bool_t strict) {
1055 Bool_t goodev = true;
1056 //
1057 if (calo2_obj && !calo2_obj->IsGood(strict))
1058 goodev = false;
1059 //
1060 if (strict) {
1061 if (trk2_obj && trk2_obj->UnpackError() != 0)
1062 goodev = false;
1063 if (tof2_obj && tof2_obj->unpackError != 0)
1064 goodev = false;
1065 if (trig_obj && trig_obj->unpackError != 0)
1066 goodev = false;
1067 if (s4_obj && s4_obj->unpackError != 0)
1068 goodev = false;
1069 if (nd_obj && nd_obj->unpackError != 0)
1070 goodev = false;
1071 if (ac_obj && (ac_obj->unpackError != 0 || ((ac_obj->status[0] >> 2) & 1) || ((ac_obj->status[1] >> 2) & 1)))
1072 goodev = false;
1073 // if(orb2_obj)
1074 }
1075 else {
1076 if (nd_obj && nd_obj->unpackError != 0)
1077 goodev = false;
1078 if (ac_obj && (ac_obj->unpackError != 0 || ((ac_obj->status[0] >> 2) & 1) || ((ac_obj->status[1] >> 2) & 1)))
1079 goodev = false;
1080 };
1081 return (goodev);
1082 }
1083 ;
1084
1085 void PamLevel2::SkipRunInfoUpdate(){
1086 printf("\n\n ******** WARNING ******** \n Skip DB connections, DO NOT USE PamLevel2::GetRunInfo() method! \n\n");
1087 fUpdateRunInfo = false;
1088 this->SetSELLI(2);
1089 printf(" ===============> W A R N I N G <================ \n");
1090 printf(" in case PamLevel2::CreateCloneTrees() will be called \n");
1091 printf(" it will be reverted to PadmeAmidala level2 structure , i.e. NO SELECTIONLIST WILL BE CREATED IN THE NEW LEVEL2 FILE! \n\n");
1092 if ( run_tree_clone ){
1093 printf(" ===============> W A R N I N G <================ \n");
1094 printf(" PamLevel2::SkipRunIndoUpdate or PamLevel2::NoDBconnections() has been called together with PamLevel2::CreateCloneTrees() \n");
1095 printf(" TO AVOID CRASHES call PamLevel2::CreateCloneTrees() after PamLevel2::SkipRunIndoUpdate or PamLevel2::NoDBconnections() \n");
1096 };
1097 }
1098
1099 void PamLevel2::SetMaxShift(Int_t sh){
1100 if ( sh >= 0 ){
1101 printf("PamLevel2::SetMaxShift(Int_t) --WARNING-- the default is optimized by checking the level2 file\n it is strongly suggested to let PamLevel2 choose the max shift!\n");
1102 maxshift = sh;
1103 } else {
1104 ULong64_t nev = GetEntries();
1105 ULong64_t runnev = 0ULL;
1106 for (Int_t r=0; r< run_tree->GetEntries();r++){
1107 run_tree->GetEntry(r);//update runinfo
1108 runnev += GetRunInfo()->NEVENTS;
1109 }
1110 maxshift = (Int_t)(runnev-nev) + 10; // +10 just to be conservative
1111 if ( (runnev-nev) == 0 ) isSync = true;
1112 if (DBG) printf("PamLevel2::SetMaxShift(Int_t_) - sh negative %i - nev is %lld runnnev is %lld so maxshift set to %i \n",sh,nev,runnev,maxshift);
1113 // printf("PamLevel2::SetMaxShift(Int_t_) - sh negative %i - nev is %lld runnnev is %lld so maxshift set to %i \n",sh,nev,runnev,maxshift); // TOGLITOGLI
1114 }
1115 }
1116
1117 //--------------------------------------
1118 //
1119 //
1120 //--------------------------------------
1121 void *PamLevel2::GetPointerTo(const char* c) {
1122
1123 TString objname = c;
1124
1125 if (!objname.CompareTo("TrkLevel1")) {
1126 if (!trk1_obj) {
1127 trk1_obj = new TrkLevel1();
1128 trk1_obj->Set();
1129 }
1130 return &trk1_obj;
1131 };
1132 if (!objname.CompareTo("TrkLevel2")) {
1133 if (!trk2_obj) {
1134 trk2_obj = new TrkLevel2();
1135 trk2_obj->Set();
1136 }
1137 return &trk2_obj;
1138 };
1139 if (!objname.CompareTo("TrkHough")) {
1140 if (!trkh_obj) {
1141 trkh_obj = new TrkHough();
1142 trkh_obj->Set();
1143 }
1144 return &trkh_obj;
1145 };
1146 if (!objname.CompareTo("CaloLevel1")) {
1147 if (!calo1_obj)
1148 calo1_obj = new CaloLevel1();
1149 return &calo1_obj;
1150 };
1151 if (!objname.CompareTo("CaloLevel2")) {
1152 if (!calo2_obj) {
1153 calo2_obj = new CaloLevel2();
1154 calo2_obj->Set();
1155 };
1156 return &calo2_obj;
1157 };
1158 if (!objname.CompareTo("ToFLevel2")) {
1159 if (!tof2_obj) {
1160 tof2_obj = new ToFLevel2();
1161 tof2_obj->Set();
1162 }
1163 return &tof2_obj;
1164 };
1165 if (!objname.CompareTo("TrigLevel2")) {
1166 if (!trig_obj)
1167 trig_obj = new TrigLevel2();
1168 return &trig_obj;
1169 };
1170 if (!objname.CompareTo("S4Level2")) {
1171 if (!s4_obj)
1172 s4_obj = new S4Level2();
1173 return &s4_obj;
1174 };
1175 if (!objname.CompareTo("NDLevel2")) {
1176 if (!nd_obj)
1177 nd_obj = new NDLevel2();
1178 return &nd_obj;
1179 };
1180 if (!objname.CompareTo("AcLevel2")) {
1181 if (!ac_obj)
1182 ac_obj = new AcLevel2();
1183 return &ac_obj;
1184 };
1185 if (!objname.CompareTo("OrbitalInfo")) {
1186 if (!orb2_obj) {
1187 orb2_obj = new OrbitalInfo();
1188 orb2_obj->Set();
1189 }
1190 return &orb2_obj;
1191 };
1192 // if(!objname.CompareTo("OrbitalInfo")){
1193 // if(!orb2_obj) orb2_obj = new OrbitalInfo();
1194 // return &orb2_obj;
1195 // };
1196 if (!objname.CompareTo("GPamela")) {
1197 if (!gp_obj)
1198 gp_obj = new GPamela();
1199 return &gp_obj;
1200 };
1201
1202 if (!objname.CompareTo("RunInfo"))
1203 return &run_obj;
1204
1205 if (!objname.CompareTo("SoftInfo"))
1206 return &soft_obj; // Emiliano
1207
1208 if (!objname.CompareTo("ProcInfo")){
1209 if (!proc_obj)
1210 proc_obj = new ProcInfo();
1211 return &proc_obj; // Emiliano
1212 }
1213
1214 return NULL;
1215 }
1216 ;
1217 //--------------------------------------
1218 //
1219 //
1220 //--------------------------------------
1221 /**
1222 * Retrieves the calorimeter track matching the seqno-th tracker stored track.
1223 * (If seqno = -1 retrieves the self-trigger calorimeter track)
1224 */
1225 CaloTrkVar *PamLevel2::GetCaloStoredTrack(int seqno) {
1226
1227 if (!calo2_obj)
1228 return 0;
1229
1230 if (calo2_obj->CaloLevel2::ntrk() == 0) {
1231 cout << "PamLevel2::GetCaloStoredTrack(int) : requested tracker SeqNo " << seqno
1232 << " but no Calorimeter tracks are stored" << endl;
1233 return NULL;
1234 };
1235
1236 CaloTrkVar *c = 0;
1237 Int_t it_calo = 0;
1238
1239 do {
1240 c = calo2_obj->CaloLevel2::GetCaloTrkVar(it_calo);
1241 it_calo++;
1242 } while (c && seqno != c->trkseqno && it_calo < calo2_obj->CaloLevel2::ntrk());
1243
1244 if (!c || seqno != c->trkseqno) {
1245 c = 0;
1246 if (seqno != -1)
1247 cout << "PamLevel2::GetCaloStoredTrack(int) : requested tracker SeqNo " << seqno
1248 << " does not match Calorimeter stored tracks" << endl;
1249 };
1250 return c;
1251
1252 }
1253 ;
1254 //--------------------------------------
1255 //
1256 //
1257 //--------------------------------------
1258 /**
1259 * Retrieves the ToF track matching the seqno-th tracker stored track.
1260 * (If seqno = -1 retrieves the tracker-independent tof track)
1261 */
1262 ToFTrkVar *PamLevel2::GetToFStoredTrack(int seqno) {
1263
1264 if (!tof2_obj)
1265 return 0;
1266
1267 if (tof2_obj->ToFLevel2::ntrk() == 0) {
1268 cout << "PamLevel2::GetToFStoredTrack(int) : requested tracker SeqNo " << seqno << " but no ToF tracks are stored"
1269 << endl;
1270 return NULL;
1271 };
1272
1273 ToFTrkVar *c = 0;
1274 Int_t it_tof = 0;
1275
1276 do {
1277 c = tof2_obj->ToFLevel2::GetToFTrkVar(it_tof);
1278 it_tof++;
1279 } while (c && seqno != c->trkseqno && it_tof < tof2_obj->ToFLevel2::ntrk());
1280
1281 if (!c || seqno != c->trkseqno) {
1282 c = 0;
1283 if (seqno != -1)
1284 cout << "PamLevel2::GetToFStoredTrack(int) : requested tracker SeqNo " << seqno
1285 << " does not match ToF stored tracks" << endl;
1286 };
1287 return c;
1288
1289 }
1290 ;
1291
1292 //--------------------------------------
1293 //
1294 //
1295 //--------------------------------------
1296 /**
1297 * Retrieves the OrbitalInfo track matching the seqno-th tracker stored track.
1298 * (If seqno = -1 retrieves the tracker-independent tof related track)
1299 */
1300 OrbitalInfoTrkVar *PamLevel2::GetOrbitalInfoStoredTrack(int seqno) {
1301
1302 if (!orb2_obj)
1303 return 0;
1304
1305 if (orb2_obj->OrbitalInfo::ntrk() == 0) {
1306 // // TRICK BEGIN
1307 // OrbitalInfoTrkVar *r = new OrbitalInfoTrkVar(); // TEMPORARY TRICK
1308 // Int_t nn = 0;
1309 // TClonesArray &tor = *orb2_obj->OrbitalInfoTrk;
1310 // for(Int_t nt=0; nt < tof2_obj->ToFLevel2::ntrk(); nt++){
1311 // //
1312 // ToFTrkVar *ptt = tof2_obj->ToFLevel2::GetToFTrkVar(nt);
1313 // if ( ptt->trkseqno != -1 ){
1314 // //
1315 // r->trkseqno = ptt->trkseqno;
1316 // //
1317 // r->Eij = 0;
1318 // //
1319 // r->Sij = 0;
1320 // //
1321 // r->pitch = -1000.;
1322 // //
1323 // r->cutoff = -1000.;
1324 // //
1325 // new(tor[nn]) OrbitalInfoTrkVar(*r);
1326 // nn++;
1327 // //
1328 // r->Clear();
1329 // //
1330 // };
1331 // };
1332 // delete r;
1333 // OrbitalInfoTrkVar *c = 0;
1334 // c = orb2_obj->OrbitalInfo::GetOrbitalInfoTrkVar(0);
1335 // return c;
1336 // //TRICK END
1337 cout << "PamLevel2::GetOrbitalInfoStoredTrack(int) : requested tracker SeqNo " << seqno
1338 << " but no OrbitalInfo tracks are stored" << endl;
1339 return NULL;
1340 };
1341
1342 OrbitalInfoTrkVar *c = 0;
1343 Int_t it_tof = 0;
1344
1345 do {
1346 c = orb2_obj->OrbitalInfo::GetOrbitalInfoTrkVar(it_tof);
1347 it_tof++;
1348 } while (c && seqno != c->trkseqno && it_tof < orb2_obj->OrbitalInfo::ntrk());
1349
1350 if (!c || seqno != c->trkseqno) {
1351 c = 0;
1352 if (seqno != -1)
1353 cout << "PamLevel2::GetOrbitalInfoStoredTrack(int) : requested tracker SeqNo " << seqno
1354 << " does not match OrbitalInfo stored tracks" << endl;
1355 };
1356 return c;
1357
1358 }
1359 ;
1360
1361 //--------------------------------------
1362 //
1363 //
1364 //--------------------------------------
1365 /**
1366 * Give the pamela track associated to a tracker track, retrieving related calorimeter, orbitalinfo and tof track information.
1367 */
1368 PamTrack* PamLevel2::GetPamTrackAlong(TrkTrack* t) {
1369
1370 cout << "PamLevel2::GetPamTrackAlong(TrkTrack* t) **obsolete** " << endl;
1371 cout << "(if you use it, remember to delete the PamTrack object)" << endl;
1372
1373 CaloTrkVar *c = 0;
1374 ToFTrkVar *o = 0;
1375 OrbitalInfoTrkVar *r = 0;
1376
1377 if (CAL2)
1378 c = GetCaloStoredTrack(t->GetSeqNo());
1379 if (TOF)
1380 o = GetToFStoredTrack(t->GetSeqNo());
1381 if (ORB)
1382 r = GetOrbitalInfoStoredTrack(t->GetSeqNo());
1383
1384 // if(t && c && o)track = new PamTrack(t,c,o);
1385 PamTrack *track = new PamTrack(t, c, o, r);
1386
1387 return track;
1388
1389 }
1390 ;
1391 //--------------------------------------
1392 //
1393 //
1394 //--------------------------------------
1395 /**
1396 * Retrieves the it-th stored track.
1397 * It override TrkLevel2::GetTrack(int it).
1398 * @param itrk Track number, ranging from 0 to GetNTracks().
1399 */
1400
1401 PamTrack* PamLevel2::GetStoredTrack(Int_t itrk) {
1402
1403 cout << "PamLevel2::GetStoredTrack(Int_t itrk) **to-be-updated** " << endl;
1404 cout
1405 << "for the moment, better use separately the methods: TrkLevel2::GetStoredTrack(seqno) CaloLevel2::GetCaloTrkVar(Int_t notrack) ToFLevel2::GetToFTrkVar(Int_t notrack) OrbitalInfo::GetOrbitalInfoTrkVar(Int_t notrack)"
1406 << endl;
1407 cout << "(if you use it, remember to delete the PamTrack object)" << endl;
1408 PamTrack *track = 0;
1409
1410 if (itrk >= 0 && itrk < trk2_obj->TrkLevel2::ntrk()) {
1411
1412 TrkTrack *t = trk2_obj->TrkLevel2::GetStoredTrack(itrk);
1413 track = GetPamTrackAlong(t);
1414
1415 }
1416 else {
1417 cout << "PamLevel2::GetStoredTrack(int) : tracker track SeqNo " << itrk << " does not exist (GetNTracks() = "
1418 << trk2_obj->TrkLevel2::GetNTracks() << ")" << endl;
1419 };
1420
1421 return track;
1422
1423 }
1424 //--------------------------------------
1425 //
1426
1427 /**
1428 * Sort physical (tracker) tracks. Left here as backward compatibility method.
1429 **/
1430 void PamLevel2::SortTracks(TString how) {
1431 printf(
1432 " WARNING! obsolete, use SortTracks() and SetSortingMethod(TString) instead! \n Setting sorting method to %s \n",
1433 how.Data());
1434 howtosort = how;
1435 SortTracks();
1436 }
1437 ;
1438
1439 //
1440 //--------------------------------------
1441 /**
1442 * Sort physical (tracker) tracks.
1443 * @param how String to set the sorting cryterium (es: "CAL" or "TRK+CAL+TOF" ecc...).
1444 * Sorting cryteria:
1445 * TRK: lower chi**2
1446 * CAL: lower Y spatial residual on the first calorimeter plane
1447 * TOF: bigger numebr of hit PMTs along the track, on S12 S21 S32 (where paddles are along the Y axis).
1448 * S1: (ask Emiliano)
1449 * S2: (ask Emiliano)
1450 * S3: (ask Emiliano)
1451 * GP: more GP hits
1452 * The default sorting cryterium is "TOF+CAL".
1453 *
1454 * The total number of physical tracks is always given by GetNTracks() and the it-th physical track can be retrieved by means of the methods GetTrack(int it) and GetTrack(int it, TString how).
1455 */
1456 void PamLevel2::SortTracks() {
1457
1458 //Check if the current event has already been sorted
1459 if (issorted && lastsorted == GetReadEntry()) {
1460 return;
1461 }
1462 //Reset the sort flags, just in case something will go wrong...
1463 issorted = false;
1464 lastsorted = -1;
1465
1466 TString how = howtosort;
1467
1468 // cout <<" PamLevel2::SortTracks(TString how) "<<endl;
1469 if (!trk2_obj) {
1470 cout << "void PamLevel2::SortTracks(): TrkLevel2 not loaded !!!";
1471 return;
1472 };
1473 //Save current Object count
1474 Int_t ObjectNumber = TProcessID::GetObjectCount();
1475
1476 // create TCloneArrays to store tracks and its images
1477 if (!tsorted)
1478 tsorted = new TClonesArray("PamTrack", trk2_obj->GetNTracks());
1479 tsorted->Clear("C+C");//Delete();
1480 TClonesArray &ttsorted = *tsorted;
1481
1482 if (!timage)
1483 timage = new TClonesArray("PamTrack", trk2_obj->GetNTracks());
1484 timage->Clear("C+C");//Delete();
1485 TClonesArray &ttimage = *timage;
1486
1487
1488
1489 //--------------------------------------------------
1490 // retrieve sorting method
1491 //--------------------------------------------------
1492 Bool_t use_TRK = how.Contains("TRK", TString::kIgnoreCase);
1493 Bool_t use_CAL = how.Contains("CAL", TString::kIgnoreCase);
1494 Bool_t use_TOF = how.Contains("TOF", TString::kIgnoreCase);
1495 Bool_t use_S1 = how.Contains("S1", TString::kIgnoreCase);
1496 Bool_t use_S2 = how.Contains("S2", TString::kIgnoreCase);
1497 Bool_t use_S3 = how.Contains("S3", TString::kIgnoreCase);
1498 Bool_t use_GP = how.Contains("GP", TString::kIgnoreCase);
1499
1500 if (use_TOF) {
1501 use_S1 = true;
1502 use_S2 = true;
1503 use_S3 = true;
1504 };
1505 if (!CAL2 && use_CAL)
1506 use_CAL = false;
1507 if (!TOF) {
1508 use_TOF = false;
1509 use_S1 = false;
1510 use_S2 = false;
1511 use_S3 = false;
1512 }
1513 if (!GP) {
1514 use_GP = false;
1515 }
1516
1517 if (!TRK2) {
1518 cout << "SortTracks() : without tracker does not work!!! (not yet)" << endl;
1519 return;
1520 };
1521
1522 // cout << "use_CAL "<<use_CAL<<" use_TOF "<<use_TOF<<" use_TRK "<<use_TRK <<endl;
1523
1524
1525 //--------------------------------------------------
1526 // loop over "physical" tracks sorted by the tracker
1527 //--------------------------------------------------
1528 for (Int_t i = 0; i < trk2_obj->TrkLevel2::GetNTracks(); i++) {
1529
1530 TrkTrack *ts = 0;
1531 CaloTrkVar *cs = 0;
1532 ToFTrkVar *os = 0;
1533 OrbitalInfoTrkVar *rs = 0;
1534
1535 // get tracker tracks
1536 TrkTrack *tp = trk2_obj->TrkLevel2::GetTrack(i); //tracker
1537 CaloTrkVar *cp = GetCaloStoredTrack(tp->GetSeqNo());
1538 ToFTrkVar *op = GetToFStoredTrack(tp->GetSeqNo());
1539 OrbitalInfoTrkVar *rp = GetOrbitalInfoStoredTrack(tp->GetSeqNo());
1540
1541 TrkTrack *ti = 0; //tracker (image)
1542 CaloTrkVar *ci = 0;
1543 ToFTrkVar *oi = 0;
1544 OrbitalInfoTrkVar *ri = 0;
1545 // cout << "trk track n. "<<i << " "<<hex<< tp <<dec<< endl;
1546 // if track has an image, check image selection
1547
1548 Int_t tp_score = 0; //main track sorted by the tracker
1549 Int_t ti_score = 0; //image track
1550 Int_t totp_score = 0; //main track sorted by the tracker
1551 Int_t toti_score = 0; //image track
1552
1553 if (tp->HasImage()) {
1554
1555 ti = trk2_obj->TrkLevel2::GetTrackImage(i); //tracker (image)
1556 ci = GetCaloStoredTrack(ti->GetSeqNo());
1557 oi = GetToFStoredTrack(ti->GetSeqNo());
1558 ri = GetOrbitalInfoStoredTrack(ti->GetSeqNo());
1559
1560 // cout << "its image "<<i << " "<<hex<< ti <<dec<< endl;
1561
1562 //assign starting scores
1563 tp_score = 0; //main track sorted by the tracker
1564 ti_score = 0; //image track
1565
1566 // -----------------------------------------------------------------------------------------
1567 // *****************************************************************************************
1568 // -----------------------------------------------------------------------------------------
1569 // calorimeter check
1570 // -----------------------------------------------------------------------------------------
1571 if (use_CAL && !calo2_obj) {
1572 cout << "void PamLevel2::SortTracks(): howtosort= " << how << " but CaloLevel2 not loaded !!!";
1573 return;
1574 };
1575 if (use_CAL && !cp && ci) {
1576 ti_score++;
1577 toti_score++;
1578 };
1579 if (use_CAL && cp && !ci) {
1580 tp_score++;
1581 totp_score++;
1582 };
1583 if (use_CAL && cp && ci && true) {
1584
1585 if (cp->npresh > ci->npresh && true) {
1586 tp_score++;
1587 totp_score++;
1588 };
1589 if (cp->npresh < ci->npresh && true) {
1590 ti_score++;
1591 toti_score++;
1592 };
1593
1594 // cout << "CALO "<<tp_score<<ti_score<<endl;
1595
1596 };
1597 // -----------------------------------------------------------------------------------------
1598 // *****************************************************************************************
1599 // -----------------------------------------------------------------------------------------
1600 // TOF check
1601 // -----------------------------------------------------------------------------------------
1602 // check the number of hit pmts along the track
1603 // on S12 S21 and S32, where paddles are parallel to Y axis
1604 if ((use_TOF || use_S1 || use_S2 || use_S3) && !tof2_obj) {
1605 cout << "void PamLevel2::SortTracks(): howtosort= " << how << " but ToFLevel2 not loaded !!!";
1606 return;
1607 };
1608 //
1609 if ((use_TOF || use_S1 || use_S2 || use_S3) && !op && oi) {
1610 ti_score++;
1611 toti_score++;
1612 };
1613 if ((use_TOF || use_S1 || use_S2 || use_S3) && op && !oi) {
1614 tp_score++;
1615 totp_score++;
1616 };
1617 if ((use_TOF || use_S1 || use_S2 || use_S3) && op && oi) {
1618 //
1619 Float_t sen = 0.;
1620 for (Int_t ih = 0; ih < op->npmtadc; ih++) {
1621 Int_t pl = tof2_obj->GetPlaneIndex((op->pmtadc).At(ih));
1622 if (pl == 2 || pl == 3 || pl == 4 || pl == 5)
1623 sen += (op->dedx).At(ih);
1624 };
1625 for (Int_t ih = 0; ih < oi->npmtadc; ih++) {
1626 Int_t pl = tof2_obj->GetPlaneIndex((oi->pmtadc).At(ih));
1627 if (pl == 2 || pl == 3 || pl == 4 || pl == 5)
1628 sen += (oi->dedx).At(ih);
1629 };
1630 //
1631 if (sen >= sortthr && false) { // temporary disabled NUCLEI special algorithm since the new one should work for every particle (to be checked!)
1632 //printf(" IS A NUCLEUS! en = %f \n",sen);
1633 //
1634 // is a nucleus use a different algorithm
1635 //
1636 Int_t nz = 6;
1637 Float_t zin[6]; // << define TOF z-coordinates
1638 for (Int_t ip = 0; ip < nz; ip++)
1639 zin[ip] = tof2_obj->ToFLevel2::GetZTOF(tof2_obj->ToFLevel2::GetToFPlaneID(ip)); // << read ToF plane z-coordinates
1640 Trajectory *tr = new Trajectory(nz, zin);
1641 //
1642 Int_t nphit_p = 0;
1643 Int_t nphit_i = 0;
1644 Float_t enhit_p = 0.;
1645 Float_t enhit_i = 0.;
1646 //
1647 for (Int_t ih = 0; ih < op->npmtadc; ih++) {
1648 Int_t pl = tof2_obj->GetPlaneIndex((op->pmtadc).At(ih));
1649 if (pl == 1 || pl == 2 || pl == 5) {
1650 nphit_p++;
1651 enhit_p += (op->dedx).At(ih);
1652 };
1653 };
1654 //
1655 tp->DoTrack2(tr);
1656 //
1657 if (fabs(tr->y[0] - oi->ytofpos[0]) < 2.) {
1658 for (Int_t ih = 0; ih < op->npmtadc; ih++) {
1659 Int_t pl = tof2_obj->GetPlaneIndex((op->pmtadc).At(ih));
1660 if (pl == 0) {
1661 nphit_p++;
1662 enhit_p += (op->dedx).At(ih);
1663 };
1664 };
1665 };
1666 if (fabs(tr->y[3] - oi->ytofpos[1]) < 2.) {
1667 for (Int_t ih = 0; ih < op->npmtadc; ih++) {
1668 Int_t pl = tof2_obj->GetPlaneIndex((op->pmtadc).At(ih));
1669 if (pl == 3) {
1670 nphit_p++;
1671 enhit_p += (op->dedx).At(ih);
1672 };
1673 };
1674 };
1675 if (fabs(tr->y[4] - oi->ytofpos[2]) < 2.) {
1676 for (Int_t ih = 0; ih < op->npmtadc; ih++) {
1677 Int_t pl = tof2_obj->GetPlaneIndex((op->pmtadc).At(ih));
1678 if (pl == 4) {
1679 nphit_p++;
1680 enhit_p += (op->dedx).At(ih);
1681 };
1682 };
1683 };
1684
1685 for (Int_t ih = 0; ih < oi->npmtadc; ih++) {
1686 Int_t pl = tof2_obj->GetPlaneIndex((oi->pmtadc).At(ih));
1687 if (pl == 1 || pl == 2 || pl == 5) {
1688 nphit_i++;
1689 enhit_i += (op->dedx).At(ih);
1690 };
1691 };
1692 //
1693 ti->DoTrack2(tr);
1694 //
1695 if (fabs(tr->y[0] - oi->ytofpos[0]) < 2.) {
1696 for (Int_t ih = 0; ih < oi->npmtadc; ih++) {
1697 Int_t pl = tof2_obj->GetPlaneIndex((oi->pmtadc).At(ih));
1698 if (pl == 0) {
1699 nphit_i++;
1700 enhit_i += (op->dedx).At(ih);
1701 };
1702 };
1703 };
1704 if (fabs(tr->y[3] - oi->ytofpos[1]) < 2.) {
1705 for (Int_t ih = 0; ih < oi->npmtadc; ih++) {
1706 Int_t pl = tof2_obj->GetPlaneIndex((oi->pmtadc).At(ih));
1707 if (pl == 3) {
1708 nphit_i++;
1709 enhit_i += (op->dedx).At(ih);
1710 };
1711 };
1712 };
1713 if (fabs(tr->y[4] - oi->ytofpos[2]) < 2.) {
1714 for (Int_t ih = 0; ih < oi->npmtadc; ih++) {
1715 Int_t pl = tof2_obj->GetPlaneIndex((oi->pmtadc).At(ih));
1716 if (pl == 4) {
1717 nphit_i++;
1718 enhit_i += (op->dedx).At(ih);
1719 };
1720 };
1721 };
1722
1723 if ((use_TOF || use_S1 || use_S2 || use_S3) && (nphit_p + nphit_i) != 0 && true) {
1724
1725 // printf(" seqno %i nphit_p %i nphit_i %i enhit_p %f enhit_i %f \n",trk2_obj->TrkLevel2::GetSeqNo(i),nphit_p,nphit_i,enhit_p,enhit_i);
1726 // printf(" score p %i score i %i \n",tp_score,ti_score);
1727 // if( enhit_p > enhit_i ) tp_score++;
1728 // if( nphit_p >= nphit_i && enhit_p > enhit_i ) tp_score++;
1729 if (nphit_p > nphit_i)
1730 tp_score++;
1731 if (nphit_p < nphit_i)
1732 ti_score++;
1733 if (nphit_p == nphit_i) {
1734 if (enhit_p > enhit_i)
1735 tp_score++;
1736 else
1737 ti_score++;
1738 };
1739 // printf(" dopo score p %i score i %i \n",tp_score,ti_score);
1740 };
1741 delete tr;
1742 //
1743 }
1744 else {
1745 // -------------
1746 // NOT a NUCLEUS
1747 // -------------
1748 //printf(" NOT a NUCLEUS! en = %f \n",sen);
1749
1750 Int_t nphit_p = 0;
1751 Int_t nphit_i = 0;
1752
1753 /* cout << "track: npmtadc "<< op->npmtadc << endl;
1754 cout << "track: npmttdc "<< op->npmttdc << endl;
1755 cout << "image: npmtadc "<< oi->npmtadc << endl;
1756 cout << "image: npmttdc "<< oi->npmttdc << endl;*/
1757
1758 // for (Int_t ih=0; ih < op->npmtadc; ih++){
1759 // Int_t pl = tof2_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
1760 // if(pl == 1 || pl == 2 || pl == 5)nphit_p++;
1761 // };
1762
1763 // for (Int_t ih=0; ih < oi->npmtadc; ih++){
1764 // Int_t pl = tof2_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
1765 // if(pl == 1 || pl == 2 || pl == 5)nphit_i++;
1766 // };
1767 // --- modified to count tdc signals (more efficient?)
1768 // --- and to implement check on tdcflag
1769 for (Int_t ih = 0; ih < op->npmttdc; ih++) {
1770 Int_t pl = tof2_obj->GetPlaneIndex((op->pmttdc).At(ih));
1771 // if( (op->tdcflag).At(ih)==0 && (pl == 1 || pl == 2 || pl == 5) )nphit_p++;
1772 if ((use_S1 && (pl == 0 || pl == 1)) || (use_S2 && (pl == 2 || pl == 3))
1773 || (use_S3 && (pl == 4 || pl == 5))) {
1774 if ((op->tdcflag).At(ih) == 0)
1775 nphit_p++;
1776 };
1777 };
1778
1779 for (Int_t ih = 0; ih < oi->npmttdc; ih++) {
1780 Int_t pl = tof2_obj->GetPlaneIndex((oi->pmttdc).At(ih));
1781 // if( (oi->tdcflag).At(ih)==0 && (pl == 1 || pl == 2 || pl == 5) )nphit_i++;
1782 if ((use_S1 && (pl == 0 || pl == 1)) || (use_S2 && (pl == 2 || pl == 3))
1783 || (use_S3 && (pl == 4 || pl == 5))) {
1784 if ((oi->tdcflag).At(ih) == 0)
1785 nphit_i++;
1786 };
1787 };
1788
1789 if ((nphit_p + nphit_i) != 0 && true) {
1790
1791 if (nphit_p != nphit_i) {
1792 totp_score += nphit_p;
1793 toti_score += nphit_i;
1794 tp_score += nphit_p;
1795 ti_score += nphit_i;
1796 };
1797 // if ( nphit_p > nphit_i) tp_score+=nphit_p;
1798 // else if( nphit_p < nphit_i) ti_score+=nphit_i;
1799 // else ;//niente
1800 };
1801 };
1802 // cout << "TOF "<<tp_score<<ti_score<<endl;
1803 };
1804
1805 // if(tp_score == ti_score) use_TRK = true;
1806
1807
1808 // -----------------------------------------------------------------------------------------
1809 // *****************************************************************************************
1810 // -----------------------------------------------------------------------------------------
1811 // tracker check
1812 // -----------------------------------------------------------------------------------------
1813 // chi**2 difference is not always large enough to distinguish among
1814 // the real track and its image.
1815 // Tracker check will be applied always when calorimeter and tof information is ambiguous.
1816 // *** MODIFIED ON AUGUST 2007 ***
1817 if (use_TRK) {
1818 // if( tp->chi2 > 0 && tp->chi2 < ti->chi2 ) tp_score++ ;
1819 // else if( ti->chi2 > 0 && ti->chi2 < tp->chi2 ) ti_score++ ;
1820
1821 // CHECK 1 : number of points along X
1822 if (tp->GetNX() >= ti->GetNX()) {
1823 tp_score++;
1824 totp_score++;
1825 };
1826 if (tp->GetNX() <= ti->GetNX()) {
1827 ti_score++;
1828 toti_score++;
1829 };
1830 // CHECK 2 : number of points along Y
1831 if (tp->GetNY() >= ti->GetNY()) {
1832 tp_score++;
1833 totp_score++;
1834 };
1835 if (tp->GetNY() <= ti->GetNY()) {
1836 ti_score++;
1837 toti_score++;
1838 };
1839
1840 // cout << "TRK "<<tp_score<<ti_score<<endl;
1841 };
1842
1843 // -----------------------------------------------------------------------------------------
1844 // *****************************************************************************************
1845 // -----------------------------------------------------------------------------------------
1846 // GPamela check
1847 // -----------------------------------------------------------------------------------------
1848
1849 //---------------------------------------------------
1850 // count the number of GP hits
1851 //---------------------------------------------------
1852 if (use_GP) {
1853 int ngphits_p = 0;
1854 int ngphits_i = 0;
1855 float toll = 0.02; //200 micron
1856 for (int ih = 0; ih < GetGPamela()->Nthspe; ih++) {
1857 int ip = (Int_t) GetGPamela()->Itrpb[ih] - 1;
1858 if (tp && tp->YGood(ip) && fabs(tp->ym[ip] - GetGPamela()->Yavspe[ih]) < toll && true)
1859 ngphits_p++;
1860 if (ti && ti->YGood(ip) && fabs(ti->ym[ip] - GetGPamela()->Yavspe[ih]) < toll && true)
1861 ngphits_i++;
1862 }
1863 if (ngphits_p > ngphits_i && true) {
1864 tp_score++;
1865 totp_score++;
1866 }
1867 if (ngphits_p < ngphits_i && true) {
1868 ti_score++;
1869 toti_score++;
1870 }
1871 }
1872
1873 // *-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
1874 // the winner is....
1875 // *-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
1876 if (tp_score > ti_score) {
1877
1878 }
1879 else if (tp_score < ti_score) {
1880
1881 ts = ti;//its image!!
1882 cs = ci;
1883 os = oi;
1884 rs = ri;
1885 Int_t totis = toti_score;
1886
1887 ti = tp;//its image!!
1888 ci = cp;
1889 oi = op;
1890 ri = rp;
1891
1892 tp = ts;//its image!!
1893 cp = cs;
1894 op = os;
1895 rp = rs;
1896
1897 toti_score = totp_score;
1898 totp_score = totis;
1899
1900 }
1901 else {
1902
1903 // cout << "Warning - track image ambiguity not solved" << endl;
1904
1905 };
1906
1907 }
1908 else {
1909 totp_score = 1;
1910 toti_score = 0;
1911
1912 // ts = tp;
1913 // cs = cp;
1914 // os = op;
1915 };
1916
1917 // cout <<" SortTracks() "<<i<<" -- "<<ts<<endl;
1918 // sorted_tracks->Add(ts);//save the track in the sorted array
1919 // sorted_tracks.Add(ts);//save the track in the sorted array
1920 // sorted_tracks.Add(tp);//save the track in the sorted array
1921 // cout << "SortTracks:: sorted_tracks->Add(it) "<<i<<" "<<ts<<endl;
1922 // cout<<"o "<<tp<<endl;
1923 // cout<<"o "<<cp<<endl;
1924 // cout<<"o "<<op<<endl;
1925
1926 new (ttsorted[i]) PamTrack(tp, cp, op, rp);
1927 new (ttimage[i]) PamTrack(ti, ci, oi, ri);
1928
1929 ((PamTrack*) (ttsorted[i]))->SetPScore(totp_score);
1930 ((PamTrack*) (ttsorted[i]))->SetIScore(toti_score);
1931 ((PamTrack*) (ttimage[i]))->SetPScore(totp_score);
1932 ((PamTrack*) (ttimage[i]))->SetIScore(toti_score);
1933 };
1934
1935 if (tsorted->GetEntries() != trk2_obj->GetNTracks()) {
1936 cout << "void PamLevel2::SortTracks(): tsorted->GetEntries() " << tsorted->GetEntries()
1937 << " != trk2_obj->GetNTracks() = " << trk2_obj->GetNTracks() << endl;
1938 tsorted->Delete();
1939 tsorted = 0;
1940 timage->Delete();
1941 timage = 0;
1942 }
1943
1944 //Restore Object count
1945 //To save space in the table keeping track of all referenced objects
1946 //We reset the object count to what it was at the beginning of the event.
1947 TProcessID::SetObjectCount(ObjectNumber);
1948
1949 //Everything went fine so the current event can be tagged as sorted
1950 issorted = true;
1951 lastsorted = GetReadEntry();
1952
1953 }
1954 ;
1955 //
1956 //--------------------------------------
1957 /**
1958 * Sort physical (tracker) tracks.
1959 * @param how String to set the sorting cryterium (es: "CAL" or "TRK+CAL+TOF" ecc...).
1960 * Sorting cryteria:
1961 * TRK: lower chi**2
1962 * CAL: lower Y spatial residual on the first calorimeter plane
1963 * TOF: bigger numebr of hit PMTs along the track, on S12 S21 S32 (where paddles are along the Y axis).
1964 * S1: (ask Emiliano)
1965 * S2: (ask Emiliano)
1966 * S3: (ask Emiliano)
1967 * GP: more GP hits
1968 * The default sorting cryterium is "TOF+CAL".
1969 *
1970 * The total number of physical tracks is always given by GetNTracks() and the it-th physical track can be retrieved by means of the methods GetTrack(int it) and GetTrack(int it, TString how).
1971 *
1972 * New version, with handling of extended tracks and nuclei tracks.
1973 */
1974 void PamLevel2::SortTracksNew() {
1975
1976 //--------------------------------------------------
1977 //Check if the current event has already been sorted
1978 //--------------------------------------------------
1979 if (issorted_new && lastsorted_new == GetReadEntry()) {
1980 return; //already done for this event
1981 }
1982
1983
1984 // cout << "SORT" << endl;
1985
1986 //Reset the sort flags, just in case something will go wrong...
1987 issorted_new = false;
1988 lastsorted_new = -1;
1989
1990 //--------------------------------------------------
1991 // set input variables
1992 //--------------------------------------------------
1993
1994 TString how = howtosort;
1995 //Save current Object count
1996 Int_t ObjectNumber = TProcessID::GetObjectCount();
1997
1998
1999
2000
2001 TrkLevel2 * trk2 ;
2002 CaloLevel2 * calo2;
2003 ToFLevel2 * tof2 ;
2004 OrbitalInfo * orb2 ;
2005
2006 TClonesArray * trkext ;
2007 TClonesArray * caloext;
2008 TClonesArray * tofext ;
2009 TClonesArray * orbext ;
2010
2011
2012 // cout << "trk2_obj" << trk2_obj<< endl;
2013 // cout << "trk2_nuc_obj" << trk2_nuc_obj<< endl;
2014 // cout << " trk_ext_obj "<<trk_ext_obj<< endl;
2015 // cout << " trk_ext_nuc_obj "<<trk_ext_nuc_obj<< endl;
2016
2017 //-----------------------------------------------------------
2018 // create/reset TCloneArrays to store tracks and their images
2019 //-----------------------------------------------------------
2020
2021 // cout << " PamLevel2::SortTracksNew() --- Clear TClonesArray objects"<<endl;
2022
2023 // main tracks from standard alg
2024 // if (!tsorted)
2025 // tsorted = new TClonesArray("PamTrack", trk2_obj->GetNTracks());
2026 // tsorted->Clear("C+C");//Delete();
2027 // // track images from standard alg
2028 // if (!timage)
2029 // timage = new TClonesArray("PamTrack", trk2_obj->GetNTracks());
2030 // timage->Clear("C+C");//Delete();
2031 // // tracks from extended algorythm
2032 // if(EXT && !text)
2033 // text = new TClonesArray("PamTrack",trk_ext_obj->GetEntries());
2034 // if(text)text->Clear("C+C");//Delete();
2035
2036 if(tsorted)delete tsorted;
2037 if(timage) delete timage;
2038 if(text) delete text;
2039 tsorted = new TClonesArray("PamTrack", trk2_obj->GetNTracks());
2040 timage = new TClonesArray("PamTrack", trk2_obj->GetNTracks());
2041 text = new TClonesArray("PamTrack",trk_ext_obj->GetEntries());
2042
2043 //-----------------------------------------------------------
2044 // create/reset TCloneArrays to store tracks and their images
2045 //-----------------------------------------------------------
2046 if(NUC){
2047
2048
2049 if(tsorted_nuc)delete tsorted_nuc;
2050 if(timage_nuc) delete timage_nuc;
2051 if(text_nuc) delete text_nuc;
2052 tsorted_nuc = new TClonesArray("PamTrack", trk2_nuc_obj->GetNTracks());
2053 timage_nuc = new TClonesArray("PamTrack", trk2_nuc_obj->GetNTracks());
2054 text_nuc = new TClonesArray("PamTrack",trk_ext_nuc_obj->GetEntries());
2055
2056 // main tracks from standard alg
2057 // if (!tsorted_nuc)
2058 // tsorted_nuc = new TClonesArray("PamTrack", trk2_nuc_obj->GetNTracks());
2059 // tsorted_nuc->Clear("C+C");//Delete();
2060 // // track images from standard alg
2061 // if (!timage_nuc)
2062 // timage_nuc = new TClonesArray("PamTrack", trk2_nuc_obj->GetNTracks());
2063 // timage_nuc->Clear("C+C");//Delete();
2064 // // tracks from extended algorythm
2065 // if(EXT && !text_nuc)
2066 // text_nuc = new TClonesArray("PamTrack",trk_ext_nuc_obj->GetEntries());
2067 // if(text_nuc)text_nuc->Clear("C+C");//Delete();
2068
2069 }
2070 //--------------------------------------------------
2071 // retrieve sorting method
2072 //--------------------------------------------------
2073 Bool_t use_TRK = how.Contains("TRK", TString::kIgnoreCase);
2074 Bool_t use_CAL = how.Contains("CAL", TString::kIgnoreCase);
2075 Bool_t use_TOF = how.Contains("TOF", TString::kIgnoreCase);
2076 Bool_t use_S1 = how.Contains("S1", TString::kIgnoreCase);
2077 Bool_t use_S2 = how.Contains("S2", TString::kIgnoreCase);
2078 Bool_t use_S3 = how.Contains("S3", TString::kIgnoreCase);
2079 Bool_t use_GP = how.Contains("GP", TString::kIgnoreCase);
2080
2081 if (use_TOF) {
2082 use_S1 = true;
2083 use_S2 = true;
2084 use_S3 = true;
2085 };
2086 if (!CAL2 && use_CAL)
2087 use_CAL = false;
2088 if (!TOF) {
2089 use_TOF = false;
2090 use_S1 = false;
2091 use_S2 = false;
2092 use_S3 = false;
2093 }
2094 if (!GP) {
2095 use_GP = false;
2096 }
2097
2098 if (!TRK2) {
2099 cout << "SortTracks() : without tracker does not work!!! (not yet)" << endl;
2100 return;
2101 };
2102
2103
2104 ///////////////////////////////////////////////////////////////////////////////////
2105 //
2106 // sort tracks and fill PamTrack arrays
2107 //
2108 ///////////////////////////////////////////////////////////////////////////////////
2109 for(int doit=0; doit<2; doit++){
2110
2111 // cout << "doit "<<doit<<endl;
2112
2113
2114 if(doit == 0){
2115
2116 trk2 = trk2_obj;
2117 calo2 = calo2_obj;
2118 tof2 = tof2_obj;
2119 orb2 = orb2_obj;
2120
2121 trkext = trk_ext_obj;
2122 caloext = calo_ext_obj;
2123 tofext = tof_ext_obj;
2124 orbext = orb_ext_obj;
2125
2126
2127
2128
2129
2130 }else if (doit == 1){
2131
2132 if(!NUC)break;
2133
2134
2135 trk2 = trk2_nuc_obj;
2136 calo2 = calo2_nuc_obj;
2137 tof2 = tof2_nuc_obj;
2138 orb2 = orb2_nuc_obj;
2139
2140 trkext = trk_ext_nuc_obj;
2141 caloext = calo_ext_nuc_obj;
2142 tofext = tof_ext_nuc_obj;
2143 orbext = orb_ext_nuc_obj;
2144
2145
2146
2147
2148 }
2149
2150 // cout << "trk2" << trk2<<endl;
2151 // cout << "calo2" << calo2<<endl;
2152 // cout << "tof2" << tof2<<endl;
2153 // cout << "orb2" << orb2<<endl;
2154 // cout << "trkext" <<trkext <<endl;
2155 // cout << "tofext" << tofext<<endl;
2156 // cout << "caloext" << caloext<<endl;
2157 // cout << "orbext" << orbext<<endl;
2158
2159 TClonesArray &ttext = (doit==0 ? *text : *text_nuc);
2160 TClonesArray &ttimage = (doit==0 ? *timage : *timage_nuc);
2161 TClonesArray &ttsorted = (doit==0 ? *tsorted : *tsorted_nuc);
2162
2163 // cout << "tsorted + timage "<<doit<<endl;
2164
2165 //--------------------------------------------------
2166 // loop over "physical" tracks sorted by the tracker
2167 //--------------------------------------------------
2168 for (Int_t i = 0; i < trk2->TrkLevel2::GetNTracks(); i++) {
2169
2170 TrkTrack *ts = 0;
2171 CaloTrkVar *cs = 0;
2172 ToFTrkVar *os = 0;
2173 OrbitalInfoTrkVar *rs = 0;
2174
2175 // get tracker tracks
2176 TrkTrack *tp = trk2->GetTrack(i); //tracker
2177 CaloTrkVar *cp = calo2->GetCaloStoredTrack(tp->GetSeqNo());
2178 ToFTrkVar *op = tof2->GetToFStoredTrack(tp->GetSeqNo());
2179 OrbitalInfoTrkVar *rp = orb2->GetOrbitalInfoStoredTrack(tp->GetSeqNo());
2180
2181 TrkTrack *ti = 0; //tracker (image)
2182 CaloTrkVar *ci = 0;
2183 ToFTrkVar *oi = 0;
2184 OrbitalInfoTrkVar *ri = 0;
2185 // cout << "trk track n. "<<i << " "<<hex<< tp <<dec<< endl;
2186 // if track has an image, check image selection
2187
2188 Int_t tp_score = 0; //main track sorted by the tracker
2189 Int_t ti_score = 0; //image track
2190 Int_t totp_score = 0; //main track sorted by the tracker
2191 Int_t toti_score = 0; //image track
2192
2193 if (tp->HasImage()) {
2194
2195 ti = trk2->GetTrackImage(i); //tracker (image)
2196 ci = calo2->GetCaloStoredTrack(ti->GetSeqNo());
2197 oi = tof2->GetToFStoredTrack(ti->GetSeqNo());
2198 ri = orb2->GetOrbitalInfoStoredTrack(ti->GetSeqNo());
2199
2200 // cout << "its image "<<i << " "<<hex<< ti <<dec<< endl;
2201
2202 //assign starting scores
2203 tp_score = 0; //main track sorted by the tracker
2204 ti_score = 0; //image track
2205
2206 // -----------------------------------------------------------------------------------------
2207 // *****************************************************************************************
2208 // -----------------------------------------------------------------------------------------
2209 // calorimeter check
2210 // -----------------------------------------------------------------------------------------
2211 if (use_CAL && !calo2) {
2212 cout << "void PamLevel2::SortTracks(): howtosort= " << how << " but CaloLevel2 not loaded !!!";
2213 return;
2214 };
2215 if (use_CAL && !cp && ci) {
2216 ti_score++;
2217 toti_score++;
2218 };
2219 if (use_CAL && cp && !ci) {
2220 tp_score++;
2221 totp_score++;
2222 };
2223 if (use_CAL && cp && ci && true) {
2224
2225 if (cp->npresh > ci->npresh && true) {
2226 tp_score++;
2227 totp_score++;
2228 };
2229 if (cp->npresh < ci->npresh && true) {
2230 ti_score++;
2231 toti_score++;
2232 };
2233
2234 // cout << "CALO "<<tp_score<<ti_score<<endl;
2235
2236 };
2237 // -----------------------------------------------------------------------------------------
2238 // *****************************************************************************************
2239 // -----------------------------------------------------------------------------------------
2240 // TOF check
2241 // -----------------------------------------------------------------------------------------
2242 // check the number of hit pmts along the track
2243 // on S12 S21 and S32, where paddles are parallel to Y axis
2244 if ((use_TOF || use_S1 || use_S2 || use_S3) && !tof2) {
2245 cout << "void PamLevel2::SortTracks(): howtosort= " << how << " but ToFLevel2 not loaded !!!";
2246 return;
2247 };
2248 //
2249 if ((use_TOF || use_S1 || use_S2 || use_S3) && !op && oi) {
2250 ti_score++;
2251 toti_score++;
2252 };
2253 if ((use_TOF || use_S1 || use_S2 || use_S3) && op && !oi) {
2254 tp_score++;
2255 totp_score++;
2256 };
2257 if ((use_TOF || use_S1 || use_S2 || use_S3) && op && oi) {
2258 //
2259 Float_t sen = 0.;
2260 for (Int_t ih = 0; ih < op->npmtadc; ih++) {
2261 Int_t pl = tof2->GetPlaneIndex((op->pmtadc).At(ih));
2262 if (pl == 2 || pl == 3 || pl == 4 || pl == 5)
2263 sen += (op->dedx).At(ih);
2264 };
2265 for (Int_t ih = 0; ih < oi->npmtadc; ih++) {
2266 Int_t pl = tof2->GetPlaneIndex((oi->pmtadc).At(ih));
2267 if (pl == 2 || pl == 3 || pl == 4 || pl == 5)
2268 sen += (oi->dedx).At(ih);
2269 };
2270 //
2271 if (sen >= sortthr && false) { // temporary disabled NUCLEI special algorithm since the new one should work for every particle (to be checked!)
2272 //printf(" IS A NUCLEUS! en = %f \n",sen);
2273 //
2274 // is a nucleus use a different algorithm
2275 //
2276 Int_t nz = 6;
2277 Float_t zin[6]; // << define TOF z-coordinates
2278 for (Int_t ip = 0; ip < nz; ip++)
2279 zin[ip] = tof2->GetZTOF(tof2->GetToFPlaneID(ip)); // << read ToF plane z-coordinates
2280 Trajectory *tr = new Trajectory(nz, zin);
2281 //
2282 Int_t nphit_p = 0;
2283 Int_t nphit_i = 0;
2284 Float_t enhit_p = 0.;
2285 Float_t enhit_i = 0.;
2286 //
2287 for (Int_t ih = 0; ih < op->npmtadc; ih++) {
2288 Int_t pl = tof2->GetPlaneIndex((op->pmtadc).At(ih));
2289 if (pl == 1 || pl == 2 || pl == 5) {
2290 nphit_p++;
2291 enhit_p += (op->dedx).At(ih);
2292 };
2293 };
2294 //
2295 tp->DoTrack2(tr);
2296 //
2297 if (fabs(tr->y[0] - oi->ytofpos[0]) < 2.) {
2298 for (Int_t ih = 0; ih < op->npmtadc; ih++) {
2299 Int_t pl = tof2->GetPlaneIndex((op->pmtadc).At(ih));
2300 if (pl == 0) {
2301 nphit_p++;
2302 enhit_p += (op->dedx).At(ih);
2303 };
2304 };
2305 };
2306 if (fabs(tr->y[3] - oi->ytofpos[1]) < 2.) {
2307 for (Int_t ih = 0; ih < op->npmtadc; ih++) {
2308 Int_t pl = tof2->GetPlaneIndex((op->pmtadc).At(ih));
2309 if (pl == 3) {
2310 nphit_p++;
2311 enhit_p += (op->dedx).At(ih);
2312 };
2313 };
2314 };
2315 if (fabs(tr->y[4] - oi->ytofpos[2]) < 2.) {
2316 for (Int_t ih = 0; ih < op->npmtadc; ih++) {
2317 Int_t pl = tof2->GetPlaneIndex((op->pmtadc).At(ih));
2318 if (pl == 4) {
2319 nphit_p++;
2320 enhit_p += (op->dedx).At(ih);
2321 };
2322 };
2323 };
2324
2325 for (Int_t ih = 0; ih < oi->npmtadc; ih++) {
2326 Int_t pl = tof2->GetPlaneIndex((oi->pmtadc).At(ih));
2327 if (pl == 1 || pl == 2 || pl == 5) {
2328 nphit_i++;
2329 enhit_i += (op->dedx).At(ih);
2330 };
2331 };
2332 //
2333 ti->DoTrack2(tr);
2334 //
2335 if (fabs(tr->y[0] - oi->ytofpos[0]) < 2.) {
2336 for (Int_t ih = 0; ih < oi->npmtadc; ih++) {
2337 Int_t pl = tof2->GetPlaneIndex((oi->pmtadc).At(ih));
2338 if (pl == 0) {
2339 nphit_i++;
2340 enhit_i += (op->dedx).At(ih);
2341 };
2342 };
2343 };
2344 if (fabs(tr->y[3] - oi->ytofpos[1]) < 2.) {
2345 for (Int_t ih = 0; ih < oi->npmtadc; ih++) {
2346 Int_t pl = tof2->GetPlaneIndex((oi->pmtadc).At(ih));
2347 if (pl == 3) {
2348 nphit_i++;
2349 enhit_i += (op->dedx).At(ih);
2350 };
2351 };
2352 };
2353 if (fabs(tr->y[4] - oi->ytofpos[2]) < 2.) {
2354 for (Int_t ih = 0; ih < oi->npmtadc; ih++) {
2355 Int_t pl = tof2->GetPlaneIndex((oi->pmtadc).At(ih));
2356 if (pl == 4) {
2357 nphit_i++;
2358 enhit_i += (op->dedx).At(ih);
2359 };
2360 };
2361 };
2362
2363 if ((use_TOF || use_S1 || use_S2 || use_S3) && (nphit_p + nphit_i) != 0 && true) {
2364
2365 // printf(" seqno %i nphit_p %i nphit_i %i enhit_p %f enhit_i %f \n",trk2_obj->TrkLevel2::GetSeqNo(i),nphit_p,nphit_i,enhit_p,enhit_i);
2366 // printf(" score p %i score i %i \n",tp_score,ti_score);
2367 // if( enhit_p > enhit_i ) tp_score++;
2368 // if( nphit_p >= nphit_i && enhit_p > enhit_i ) tp_score++;
2369 if (nphit_p > nphit_i)
2370 tp_score++;
2371 if (nphit_p < nphit_i)
2372 ti_score++;
2373 if (nphit_p == nphit_i) {
2374 if (enhit_p > enhit_i)
2375 tp_score++;
2376 else
2377 ti_score++;
2378 };
2379 // printf(" dopo score p %i score i %i \n",tp_score,ti_score);
2380 };
2381 delete tr;
2382 //
2383 }
2384 else {
2385 // -------------
2386 // NOT a NUCLEUS
2387 // -------------
2388 //printf(" NOT a NUCLEUS! en = %f \n",sen);
2389
2390 Int_t nphit_p = 0;
2391 Int_t nphit_i = 0;
2392
2393 /* cout << "track: npmtadc "<< op->npmtadc << endl;
2394 cout << "track: npmttdc "<< op->npmttdc << endl;
2395 cout << "image: npmtadc "<< oi->npmtadc << endl;
2396 cout << "image: npmttdc "<< oi->npmttdc << endl;*/
2397
2398 // for (Int_t ih=0; ih < op->npmtadc; ih++){
2399 // Int_t pl = tof2_obj->GetPlaneIndex( (op->pmtadc).At(ih) );
2400 // if(pl == 1 || pl == 2 || pl == 5)nphit_p++;
2401 // };
2402
2403 // for (Int_t ih=0; ih < oi->npmtadc; ih++){
2404 // Int_t pl = tof2_obj->GetPlaneIndex( (oi->pmtadc).At(ih) );
2405 // if(pl == 1 || pl == 2 || pl == 5)nphit_i++;
2406 // };
2407 // --- modified to count tdc signals (more efficient?)
2408 // --- and to implement check on tdcflag
2409 for (Int_t ih = 0; ih < op->npmttdc; ih++) {
2410 Int_t pl = tof2->GetPlaneIndex((op->pmttdc).At(ih));
2411 // if( (op->tdcflag).At(ih)==0 && (pl == 1 || pl == 2 || pl == 5) )nphit_p++;
2412 if ((use_S1 && (pl == 0 || pl == 1)) || (use_S2 && (pl == 2 || pl == 3))
2413 || (use_S3 && (pl == 4 || pl == 5))) {
2414 if ((op->tdcflag).At(ih) == 0)
2415 nphit_p++;
2416 };
2417 };
2418
2419 for (Int_t ih = 0; ih < oi->npmttdc; ih++) {
2420 Int_t pl = tof2->GetPlaneIndex((oi->pmttdc).At(ih));
2421 // if( (oi->tdcflag).At(ih)==0 && (pl == 1 || pl == 2 || pl == 5) )nphit_i++;
2422 if ((use_S1 && (pl == 0 || pl == 1)) || (use_S2 && (pl == 2 || pl == 3))
2423 || (use_S3 && (pl == 4 || pl == 5))) {
2424 if ((oi->tdcflag).At(ih) == 0)
2425 nphit_i++;
2426 };
2427 };
2428
2429 if ((nphit_p + nphit_i) != 0 && true) {
2430
2431 if (nphit_p != nphit_i) {
2432 totp_score += nphit_p;
2433 toti_score += nphit_i;
2434 tp_score += nphit_p;
2435 ti_score += nphit_i;
2436 };
2437 // if ( nphit_p > nphit_i) tp_score+=nphit_p;
2438 // else if( nphit_p < nphit_i) ti_score+=nphit_i;
2439 // else ;//niente
2440 };
2441 };
2442 // cout << "TOF "<<tp_score<<ti_score<<endl;
2443 };
2444
2445 // if(tp_score == ti_score) use_TRK = true;
2446
2447
2448 // -----------------------------------------------------------------------------------------
2449 // *****************************************************************************************
2450 // -----------------------------------------------------------------------------------------
2451 // tracker check
2452 // -----------------------------------------------------------------------------------------
2453 // chi**2 difference is not always large enough to distinguish among
2454 // the real track and its image.
2455 // Tracker check will be applied always when calorimeter and tof information is ambiguous.
2456 // *** MODIFIED ON AUGUST 2007 ***
2457 if (use_TRK) {
2458 // if( tp->chi2 > 0 && tp->chi2 < ti->chi2 ) tp_score++ ;
2459 // else if( ti->chi2 > 0 && ti->chi2 < tp->chi2 ) ti_score++ ;
2460
2461 // CHECK 1 : number of points along X
2462 if (tp->GetNX() >= ti->GetNX()) {
2463 tp_score++;
2464 totp_score++;
2465 };
2466 if (tp->GetNX() <= ti->GetNX()) {
2467 ti_score++;
2468 toti_score++;
2469 };
2470 // CHECK 2 : number of points along Y
2471 if (tp->GetNY() >= ti->GetNY()) {
2472 tp_score++;
2473 totp_score++;
2474 };
2475 if (tp->GetNY() <= ti->GetNY()) {
2476 ti_score++;
2477 toti_score++;
2478 };
2479
2480 // cout << "TRK "<<tp_score<<ti_score<<endl;
2481 };
2482
2483 // -----------------------------------------------------------------------------------------
2484 // *****************************************************************************************
2485 // -----------------------------------------------------------------------------------------
2486 // GPamela check
2487 // -----------------------------------------------------------------------------------------
2488
2489 //---------------------------------------------------
2490 // count the number of GP hits
2491 //---------------------------------------------------
2492 if (use_GP) {
2493 int ngphits_p = 0;
2494 int ngphits_i = 0;
2495 float toll = 0.02; //200 micron
2496 for (int ih = 0; ih < GetGPamela()->Nthspe; ih++) {
2497 int ip = (Int_t) GetGPamela()->Itrpb[ih] - 1;
2498 if (tp && tp->YGood(ip) && fabs(tp->ym[ip] - GetGPamela()->Yavspe[ih]) < toll && true)
2499 ngphits_p++;
2500 if (ti && ti->YGood(ip) && fabs(ti->ym[ip] - GetGPamela()->Yavspe[ih]) < toll && true)
2501 ngphits_i++;
2502 }
2503 if (ngphits_p > ngphits_i && true) {
2504 tp_score++;
2505 totp_score++;
2506 }
2507 if (ngphits_p < ngphits_i && true) {
2508 ti_score++;
2509 toti_score++;
2510 }
2511 }
2512
2513 // *-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
2514 // the winner is....
2515 // *-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
2516 if (tp_score > ti_score) {
2517
2518 }
2519 else if (tp_score < ti_score) {
2520
2521 ts = ti;//its image!!
2522 cs = ci;
2523 os = oi;
2524 rs = ri;
2525 Int_t totis = toti_score;
2526
2527 ti = tp;//its image!!
2528 ci = cp;
2529 oi = op;
2530 ri = rp;
2531
2532 tp = ts;//its image!!
2533 cp = cs;
2534 op = os;
2535 rp = rs;
2536
2537 toti_score = totp_score;
2538 totp_score = totis;
2539
2540 }
2541 else {
2542
2543 // cout << "Warning - track image ambiguity not solved" << endl;
2544
2545 };
2546
2547 }
2548 else {
2549 totp_score = 1;
2550 toti_score = 0;
2551
2552 // ts = tp;
2553 // cs = cp;
2554 // os = op;
2555 };
2556
2557 // cout <<" SortTracks() "<<i<<" -- "<<ts<<endl;
2558 // sorted_tracks->Add(ts);//save the track in the sorted array
2559 // sorted_tracks.Add(ts);//save the track in the sorted array
2560 // sorted_tracks.Add(tp);//save the track in the sorted array
2561 // cout << "SortTracks:: sorted_tracks->Add(it) "<<i<<" "<<ts<<endl;
2562 // cout<<"o "<<tp<<endl;
2563 // cout<<"o "<<cp<<endl;
2564 // cout<<"o "<<op<<endl;
2565
2566 new (ttsorted[i]) PamTrack(tp, cp, op, rp);
2567 new (ttimage[i]) PamTrack(ti, ci, oi, ri);
2568
2569 ((PamTrack*) (ttsorted[i]))->SetPScore(totp_score);
2570 ((PamTrack*) (ttsorted[i]))->SetIScore(toti_score);
2571 ((PamTrack*) (ttimage[i]))->SetPScore(totp_score);
2572 ((PamTrack*) (ttimage[i]))->SetIScore(toti_score);
2573 };
2574
2575
2576
2577
2578
2579 // if (tsorted->GetEntries() != trk2->GetNTracks()) {
2580 // cout << "void PamLevel2::SortTracks(): tsorted->GetEntries() " << tsorted->GetEntries()
2581 // << " != trk2->GetNTracks() = " << trk2->GetNTracks() << endl;
2582 // tsorted->Delete();
2583 // tsorted = 0;
2584 // timage->Delete();
2585 // timage = 0;
2586 // }
2587
2588
2589 // cout << "text "<<doit<<endl;
2590
2591
2592 //--------------------------------------------------
2593 // fill array with extended tracks (this is easy...)
2594 //--------------------------------------------------
2595 if(EXT){
2596 for(int it=0; it<trkext->GetEntries(); it++){
2597
2598 new (ttext[it]) PamTrack((ExtTrack*)(*trkext)[it], (CaloTrkVar*)(*caloext)[it], (ToFTrkVar*)(*tofext)[it], (OrbitalInfoTrkVar*)(*orbext)[it]);
2599
2600 }
2601 }
2602
2603
2604
2605 };
2606
2607
2608 //Restore Object count
2609 //To save space in the table keeping track of all referenced objects
2610 //We reset the object count to what it was at the beginning of the event.
2611 TProcessID::SetObjectCount(ObjectNumber);
2612
2613 //Everything went fine so the current event can be tagged as sorted
2614 issorted_new = true;
2615 lastsorted_new = GetReadEntry();
2616
2617 // cout << " tsorted "<< tsorted << " "<<(tsorted ? tsorted->GetEntries() : 0)<<endl;
2618 // cout << " timage "<< timage << " "<<(timage ? timage->GetEntries() : 0)<<endl;
2619 // cout << " text "<< text << " "<<(text ? text->GetEntries() : 0)<<endl;
2620
2621
2622 };
2623 //--------------------------------------
2624 //
2625 //
2626 //--------------------------------------
2627 /**
2628 * This method overrides TrkLevel2::GetTracks(), where sorting is done by decreasing number of fit points and increasing chi^2.
2629 * PamLevel2::GetTracks() keeps the same track order given by TrkLevel2::GetTracks(), but checks image selection by using calorimeter and ToF tracking information.
2630 */
2631
2632 // TRefArray *PamLevel2::GetTracks(){
2633
2634 // // *-*-*-*-*-*-*-*-*-*-*-*-*
2635 // SortTracks("+CAL+TOF");
2636 // // *-*-*-*-*-*-*-*-*-*-*-*-*
2637
2638 // // return sorted_tracks;
2639 // return &sorted_tracks;
2640
2641 // };
2642
2643
2644 //--------------------------------------
2645 //
2646 //
2647 //--------------------------------------
2648
2649 /**
2650 * Retrieves the it-th Pamela "physical" track.
2651 * It override TrkLevel2::GetTrack(int it).
2652 * @param it Track number, ranging from 0 to GetNTracks().
2653 */
2654 PamTrack *PamLevel2::GetTrackOld(int it) {
2655
2656 PamTrack *track = NULL;
2657
2658 // *-*-*-*-*-*-*-*-*-*-*-*-*
2659 SortTracks();
2660 // *-*-*-*-*-*-*-*-*-*-*-*-*
2661 if (!tsorted)
2662 return track;
2663 if (!tsorted->GetEntries())
2664 return track;
2665
2666
2667 if (it >= 0 && it < trk2_obj->TrkLevel2::GetNTracks()) {
2668 track = (PamTrack*)((*tsorted)[it]);
2669 }
2670 else {
2671 cout << "PamLevel2::GetTrackOld(int) : tracker track SeqNo " << it << " does not exist (GetNTracks() = "
2672 << trk2_obj->TrkLevel2::GetNTracks() << ")" << endl;
2673 };
2674
2675 return track;
2676
2677 };
2678
2679 //PamTrack *PamLevel2::GetTrack(int it) { return GetTrack(it,trkAlg); };
2680
2681 /**
2682 * Retrieves the it-th Pamela "physical" track.
2683 * It override TrkLevel2::GetTrack(int it).
2684 * @param it Track number, ranging from 0 to GetNTracks().
2685 * @param alg Algorythm, see SetTrakingAlgorythm(char *alg) for explanation.
2686 */
2687 PamTrack *PamLevel2::GetTrack(int it, const char* alg) {
2688
2689 TString s(alg);
2690 if(!s.CompareTo("") ||!s.CompareTo("STD") )return GetTrackOld(it); //old algorythm
2691
2692
2693 SortTracksNew();
2694 // >> fill tsorted, timage and text
2695
2696 if ( ( !s.Contains("EXTF", TString::kIgnoreCase) || !EXT )){ //not forced exteded-alg requested (or exteded-data missing)
2697
2698 if( s.Contains("NUC")){
2699 if(
2700 tsorted_nuc &&
2701 it < tsorted_nuc->GetEntries() && //enough tracks found
2702 it >= 0 && //valid index
2703 true) return (PamTrack*)((*tsorted_nuc)[it]); //ok return the track
2704 }else{
2705 if(
2706 tsorted &&
2707 it < tsorted->GetEntries() && //enough tracks found
2708 it >= 0 && //valid index
2709 true )return (PamTrack*)((*tsorted)[it]); //ok return the track
2710 }
2711
2712 }
2713
2714
2715 /////////////////////////////////////////////////////////////////////////
2716 /// if requested get track from extended algorythm output
2717 /////////////////////////////////////////////////////////////////////////
2718
2719 if(s.Contains("EXT", TString::kIgnoreCase) && EXT){//if exteded-alg requested
2720
2721 if(s.Contains("NUC")){
2722 if(
2723 text_nuc &&
2724 it < text_nuc->GetEntries() && //enough tracks found
2725 it >= 0 && //valid index
2726 true) return (PamTrack*)((*text_nuc)[it]);
2727 }else{
2728 if(
2729 text &&
2730 it < text->GetEntries() && //enough tracks found
2731 it >= 0 && //valid index
2732 true) return (PamTrack*)((*text)[it]);
2733 }
2734
2735 };
2736
2737
2738 cout << "PamTrack *PamLevel2::GetTrack("<<it<<","<<alg<<") -- wrong track number or unrecognised algorithm"<<endl;
2739
2740 return NULL;
2741
2742 }
2743 ;
2744 TClonesArray *PamLevel2::GetTracks() {
2745
2746 // *-*-*-*-*-*-*-*-*-*-*-*-*
2747 SortTracks();
2748 // *-*-*-*-*-*-*-*-*-*-*-*-*
2749
2750 return tsorted;
2751
2752 }
2753 ;
2754 Int_t PamLevel2::GetNTracks(const char* alg) {
2755
2756
2757 // cout << " trk_nuc_obj->GetEntries() "<<trk_nuc_obj->GetEntries()<<" trk2_nuc_obj->GetNTracks() "<<trk2_nuc_obj->GetNTracks()<<endl;
2758
2759 TString s(alg);
2760
2761 if(!s.CompareTo("") || !s.CompareTo("STD"))return trk2_obj->TrkLevel2::GetNTracks(); //standard algorythm
2762
2763 if(s.Contains("EXTF", TString::kIgnoreCase) && EXT){
2764 if(s.Contains("NUC", TString::kIgnoreCase) && NUC)return trk_ext_nuc_obj->GetEntries();//ok
2765 return trk_ext_obj->GetEntries();//ok
2766 }
2767 if( s.Contains("EXT", TString::kIgnoreCase) && EXT) {
2768 if(s.Contains("NUC", TString::kIgnoreCase) && NUC)
2769 return (trk2_nuc_obj->TrkLevel2::GetNTracks() ? trk2_nuc_obj->TrkLevel2::GetNTracks() : trk_ext_nuc_obj->GetEntries() );
2770 return (trk2_obj->TrkLevel2::GetNTracks() ? trk2_obj->TrkLevel2::GetNTracks() : trk_ext_obj->GetEntries() );
2771 }
2772 if(s.Contains("NUC", TString::kIgnoreCase) && NUC )
2773 return trk2_nuc_obj->TrkLevel2::GetNTracks();
2774
2775 cout << "Int_t PamLevel2::GetNTracks("<<alg<<") -- unrecognised algorithm"<<endl;
2776
2777 return 0;
2778
2779 }
2780
2781
2782 //--------------------------------------
2783 //
2784 //
2785 //--------------------------------------
2786 /**
2787 * Retrieves (if present) the image of the it-th Pamela "physical" track, sorted by the method PamLevel2::SortTracks().
2788 * @param it Track number, ranging from 0 to GetNTracks().
2789 */
2790 PamTrack *PamLevel2::GetTrackImageOld(int it) {
2791
2792 // *-*-*-*-*-*-*-*-*-*-*-*-*
2793 SortTracks();
2794 // *-*-*-*-*-*-*-*-*-*-*-*-*
2795 if (!timage)
2796 return 0;
2797 if (!timage->GetEntries())
2798 return 0;
2799
2800 PamTrack *image = 0;
2801
2802 if (it >= 0 && it < trk2_obj->TrkLevel2::GetNTracks()) {
2803 TClonesArray &t = *(tsorted);
2804 PamTrack *temp = (PamTrack*) t[it];
2805 if (temp->GetTrkTrack()->HasImage()) {
2806 TClonesArray & t = *(timage);
2807 image = (PamTrack*) t[it];
2808 }
2809 else {
2810 // cout <<"PamLevel2::GetTrackImage(int) : Track SeqNo "<<it<<" does not have image"<<endl;
2811 };
2812 }
2813 else {
2814 cout << "PamLevel2::GetTrackImageOld(int) : Tracker track SeqNo " << it << " does not exist (GetNTracks() = "
2815 << trk2_obj->TrkLevel2::GetNTracks() << ")" << endl;
2816 };
2817
2818 return image;
2819 }
2820 /**
2821 * Retrieves (if present) the image of the it-th Pamela "physical" track, sorted by the method PamLevel2::SortTracks().
2822 * @param it Track number, ranging from 0 to GetNTracks().
2823 * @param alg Algorythm, see SetTrakingAlgorythm(char *alg) for explanation.
2824 */
2825 PamTrack *PamLevel2::GetTrackImage(int it, const char* alg) {
2826
2827 TString s(alg);
2828 if(!s.CompareTo("") || !s.CompareTo("STD"))return GetTrackImageOld(it); //old algorythm
2829
2830
2831 SortTracksNew();
2832 // >> fill tsorted, timage and text
2833
2834 if ( ( !s.Contains("EXTF", TString::kIgnoreCase) || !EXT )){ //not forced exteded-alg requested (or exteded-data missing)
2835
2836 if( s.Contains("NUC")){
2837 if(
2838 tsorted_nuc &&
2839 it < tsorted_nuc->GetEntries() && //enough tracks found
2840 it >= 0 && //valid index
2841 true){
2842 TClonesArray &t = *(tsorted_nuc);
2843 PamTrack *temp = (PamTrack*) t[it];
2844 if (temp->GetTrkTrack()->HasImage()) {
2845 return (PamTrack*)((*timage_nuc)[it]); //ok return the track
2846 }else{
2847 return NULL;
2848 }
2849
2850 }
2851 }else{
2852 if(
2853 tsorted &&
2854 it < tsorted->GetEntries() && //enough tracks found
2855 it >= 0 && //valid index
2856 true ){
2857 TClonesArray &t = *(tsorted);
2858 PamTrack *temp = (PamTrack*) t[it];
2859 if (temp->GetTrkTrack()->HasImage()) {
2860 return (PamTrack*)((*timage)[it]); //ok return the track
2861 }else{
2862 return NULL;
2863 }
2864 }
2865 }
2866
2867 }
2868
2869 cout << "PamTrack *PamLevel2::GetTrackImage("<<it<<","<<alg<<") -- wrong track number or unrecognised algorithm"<<endl;
2870
2871 return NULL;
2872
2873 }
2874 ;
2875
2876 //--------------------------------------
2877 //
2878 //
2879 //--------------------------------------
2880 /**
2881 * Get the Pamela detector trees in a single file and make them friends.
2882 * @param f TFile pointer
2883 * @param detlist String to select trees to be included
2884 * @return Pointer to a TTree
2885 */
2886 TTree *PamLevel2::GetPamTree(TFile *f, TString detlist) {
2887
2888 if (detlist.Contains("+AUTO", TString::kIgnoreCase)) {
2889 cout << "+AUTO" << endl;
2890 GetWhichTrees(f);
2891 };
2892 SetWhichTrees(detlist);
2893
2894 if (pam_tree) {
2895 printf("WARNING: TTree *PamLevel2::GetPamTree(TFile *fl, TString detlist) -- pam_tree already exists!\n ");
2896 return pam_tree;
2897 };
2898 //
2899
2900 cout << "TTree *PamLevel2::GetPamTree(TFile *f, TString detlist ) -- obsolte " << endl;
2901
2902 // SetWhichTrees(detlist);
2903
2904 TTree *Trout = 0;
2905
2906 TString fname = f->GetName();
2907 if (!CheckLevel2File(fname))
2908 return NULL;
2909
2910 // UInt_t *found=0;
2911
2912 cout << "GetPamTree(TFile*,TString): detector list --> ";
2913 if (TRK1)
2914 cout << "TRK1 ";
2915 if (TRK2)
2916 cout << "TRK2 ";
2917 if (TRKh)
2918 cout << "TRKH ";
2919 if (CAL1)
2920 cout << "CAL1 ";
2921 if (CAL2)
2922 cout << "CAL2 ";
2923 if (TOF)
2924 cout << "TOF ";
2925 if (TRG)
2926 cout << "TRG ";
2927 if (AC)
2928 cout << "AC ";
2929 if (ND)
2930 cout << "ND ";
2931 if (S4)
2932 cout << "S4 ";
2933 if (ORB)
2934 cout << "ORB ";
2935 if (GP)
2936 cout << "GP ";
2937 cout << endl;
2938 if (SELLI && SELLI != 2)
2939 cout << ">>> Found selection-list <<<" << endl; //EMILIANO
2940
2941 f->cd();
2942
2943 // Tracker
2944 TTree *T = (TTree*) f->Get("Tracker");
2945 if (T && (TRK2 || TRK1 || TRKh)) {
2946 if (TRK2)
2947 T->SetBranchAddress("TrkLevel2", GetPointerTo("TrkLevel2"));
2948 // else T->SetBranchStatus("TrkLevel2",0,found);
2949 if (TRK2)
2950 cout << "Tracker : set branch address TrkLevel2" << endl;
2951 if (TRK1)
2952 T->SetBranchAddress("TrkLevel1", GetPointerTo("TrkLevel1"));
2953 // else T->SetBranchStatus("TrkLevel1",0,found);
2954 if (TRK1)
2955 cout << "Tracker : set branch address TrkLevel1" << endl;
2956 if (TRKh)
2957 T->SetBranchAddress("TrkHough", GetPointerTo("TrkHough"));
2958 // else T->SetBranchStatus("TrkHough",0,found);
2959 if (TRKh)
2960 cout << "Tracker : set branch address TrkHough" << endl;
2961 if (!Trout)
2962 Trout = T;
2963 else
2964 Trout->AddFriend(T);
2965 }
2966 else {
2967 cout << "Tracker : missing tree" << endl;
2968 };
2969 // Calorimeter
2970 TTree *C = (TTree*) f->Get("Calorimeter");
2971 if (C && (CAL2 || CAL1)) {
2972 if (CAL2)
2973 C->SetBranchAddress("CaloLevel2", GetPointerTo("CaloLevel2"));
2974 // else C->SetBranchStatus("CaloLevel2",0,found);
2975 if (CAL2)
2976 cout << "Calorimeter : set branch address CaloLevel2" << endl;
2977 if (CAL1)
2978 C->SetBranchAddress("CaloLevel1", GetPointerTo("CaloLevel1"));
2979 // else C->SetBranchStatus("CaloLevel1",0,found);
2980 if (CAL1)
2981 cout << "Calorimeter : set branch address CaloLevel1" << endl;
2982 if (!Trout)
2983 Trout = C;
2984 else
2985 Trout->AddFriend(C);
2986 }
2987 else {
2988 cout << "Calorimeter : missing tree" << endl;
2989 };
2990
2991 // ToF
2992 TTree *O = (TTree*) f->Get("ToF");
2993 if (O && TOF) {
2994 O->SetBranchAddress("ToFLevel2", GetPointerTo("ToFLevel2"));
2995 cout << "ToF : set branch address ToFLevel2" << endl;
2996 if (!Trout)
2997 Trout = O;
2998 else
2999 Trout->AddFriend(O);
3000 }
3001 else {
3002 cout << "ToF : missing tree" << endl;
3003 };
3004 // Trigger
3005 TTree *R = (TTree*) f->Get("Trigger");
3006 if (R && TRG) {
3007 R->SetBranchAddress("TrigLevel2", GetPointerTo("TrigLevel2"));
3008 cout << "Trigger : set branch address TrigLevel2" << endl;
3009 if (!Trout)
3010 Trout = R;
3011 else
3012 Trout->AddFriend(R);
3013 }
3014 else {
3015 cout << "Trigger : missing tree" << endl;
3016 };
3017 // S4
3018 TTree *S = (TTree*) f->Get("S4");
3019 if (S && S4) {
3020 S->SetBranchAddress("S4Level2", GetPointerTo("S4Level2"));
3021 cout << "S4 : set branch address S4Level2" << endl;
3022 if (!Trout)
3023 Trout = S;
3024 else
3025 Trout->AddFriend(S);
3026 }
3027 else {
3028 cout << "S4 : missing tree" << endl;
3029 };
3030 // Neutron Detector
3031 TTree *N = (TTree*) f->Get("NeutronD");
3032 if (N && ND) {
3033 N->SetBranchAddress("NDLevel2", GetPointerTo("NDLevel2"));
3034 cout << "NeutronD : set branch address NDLevel2" << endl;
3035 if (!Trout)
3036 Trout = N;
3037 else
3038 Trout->AddFriend(N);
3039 }
3040 else {
3041 cout << "NeutronD : missing tree" << endl;
3042 };
3043 // Anticounters
3044 TTree *A = (TTree*) f->Get("Anticounter");
3045 if (A && AC) {
3046 A->SetBranchAddress("AcLevel2", GetPointerTo("AcLevel2"));
3047 cout << "Anticounter : set branch address AcLevel2" << endl;
3048 if (!Trout)
3049 Trout = A;
3050 else
3051 Trout->AddFriend(A);
3052 }
3053 else {
3054 cout << "Anticounter : missing tree" << endl;
3055 };
3056 // Orbital Info
3057 TTree *B = (TTree*) f->Get("OrbitalInfo");
3058 if (B && ORB) {
3059 B->SetBranchAddress("OrbitalInfo", GetPointerTo("OrbitalInfo"));
3060 cout << "OrbitalInfo : set branch address OrbitalInfo" << endl;
3061 if (!Trout)
3062 Trout = B;
3063 else
3064 Trout->AddFriend(B);
3065 }
3066 else {
3067 cout << "OrbitalInfo : missing tree" << endl;
3068 };
3069
3070 // GPamela
3071 TTree *G = (TTree*) f->Get("h20");
3072 if (G && GP) {
3073 if (!gp_obj)
3074 gp_obj = new GPamela();
3075 // ------------------------------------
3076 // ATTENZIONE!!!
3077 // non so per quale arcano motivo,
3078 // per l'albero di gpamela il branch address lo devo settare
3079 // DOPO aver fatto friend
3080 // FGRRRVZZZZUTSALKJMSLKJ!!!
3081 // ------------------------------------
3082 // gp_obj->SetBranchAddress(G); //ho dovuto fare in maniera diversa dagli altri
3083 // cout << "h20 : set branch address GPamela "<<endl;
3084 if (!Trout)
3085 Trout = G;
3086 else
3087 Trout->AddFriend(G);
3088 }
3089 else {
3090 // cout << "h20 : missing tree"<<endl;
3091 };
3092
3093 TTree *L = (TTree*) f->Get("SelectionList");
3094 if (L && SELLI == 1) {
3095 cout << " TTree *PamLevel2::GetPamTree(TFile, TString) >>> SelectionList not implemented!!!!" << endl;
3096 sel_tree = 0;
3097 }
3098 else {
3099 cout << "SelectionList : missing tree" << endl;
3100 };
3101
3102 cout << endl << " Number of entries: " << Trout->GetEntries() << endl << endl;
3103
3104 // ------------------------------------
3105 // ATTENZIONE!!!
3106 // non so per quale arcano motivo,
3107 // per l'albero di gpamela il branch address lo devo settare
3108 // DOPO aver fatto friend
3109 // FGRRRVZZZZUTSALKJMSLKJ!!!
3110 // ------------------------------------
3111 if (G && GP) {
3112 gp_obj->SetBranchAddress(Trout); //ho dovuto fare in maniera diversa dagli altri
3113 cout << "h20 : set branch address GPamela " << endl;
3114 }
3115 else {
3116 cout << "h20 : missing tree" << endl;
3117 };
3118
3119 pam_tree = (TChain*) Trout;
3120
3121 return Trout;
3122
3123 }
3124 //--------------------------------------
3125 //
3126 //
3127 //--------------------------------------
3128 /**
3129 * Get list of Level2 files.
3130 * @param ddir Level2 data directory.
3131 * @param flisttxt Name of txt file containing file list.
3132 * @return Pointer to a TList of TSystemFiles
3133 * If no input file list is given , all the Level2 files inside the directory are processed.
3134 */
3135 TList* PamLevel2::GetListOfLevel2Files(TString ddir, TString flisttxt = "") {
3136
3137 TString wdir = gSystem->WorkingDirectory();
3138
3139 if (ddir != "") {
3140 cout << "Level2 data directory : " << ddir << endl;
3141 }
3142 else {
3143 cout << "Level2 data directory not given as input: trying to evaluate from list or taking working directory "
3144 << endl;
3145 };
3146 TList *contents = new TList; // create output list
3147 contents->SetOwner();
3148
3149 // --------------------------------------
3150 // case 1 : input file/file-list provided
3151 // --------------------------------------
3152 if (flisttxt != "") {
3153
3154 // --------------------------------------
3155 // a list of files given as input
3156 // --------------------------------------
3157 if (!flisttxt.EndsWith(".root")) {
3158
3159 TString dir = gSystem->DirName(flisttxt);
3160 // cout << " List directory "<<dir<<endl;
3161 if (dir == "." || dir == "")
3162 flisttxt = gSystem->ConcatFileName(wdir.Data(), gSystem->BaseName(flisttxt));
3163 // flisttxt = gSystem->ConcatFileName(gSystem->DirName(flisttxt),gSystem->BaseName(flisttxt));
3164
3165 if (!gSystem->ChangeDirectory(ddir)) {
3166 cout << "Cannot change directory : " << ddir << endl;
3167 return 0;
3168 }
3169
3170 cout << "Input file list : " << flisttxt << endl;
3171 ifstream in;
3172 in.open(flisttxt, ios::in); //open input file list
3173 if (!in.good()) {
3174 cout << " TList* PamLevel2::GetListOfLevel2Files(TString, TString) --> ERROR opening the file " << endl;
3175 gSystem->ChangeDirectory(wdir); // back to the working directory
3176 return 0;
3177 }
3178 int line = 0;
3179 // ...........................
3180 // loop over file-list records
3181 // ...........................
3182 while (1) {
3183 TString file;
3184 in >> file;
3185 if (!in.good())
3186 break;
3187 line++;
3188 cout << "(*) " << file << endl;
3189 if (file.IsNull()) {
3190 cout << "-- list interrupted at line " << line << endl;
3191 break;
3192 }
3193 if (file.Contains("#"))
3194 file = file(0, file.First("#"));
3195 //
3196 // take only root files
3197 //
3198 if (file.EndsWith(".root")) {
3199 TString filedir;
3200 cout << ddir << endl;
3201 if ( ddir != "" ) {
3202 filedir = ddir; // take the input dir
3203 }
3204 else {
3205 gSystem->ChangeDirectory(wdir); // back to the working directory
3206 filedir = gSystem->DirName(file); // this will take the path if exist in the list otherwise it will return automatically the working dir
3207 };
3208 filedir.Append("/");
3209 // char *fullpath = gSystem->ConcatFileName(gSystem->DirName(filedir), gSystem->BaseName(file));
3210 char *fullpath = gSystem->ConcatFileName(filedir.Data(), gSystem->BaseName(file));
3211 contents->Add(new TSystemFile(fullpath, gSystem->DirName(fullpath)));// add file to the list
3212 cout << fullpath << endl;
3213 delete fullpath;
3214 }
3215 // }else{
3216 // if(file.Data()!="")cout << "File: "<<file<<" ---> missing "<< endl;
3217 // };
3218 };
3219 in.close();
3220 // --------------------------------------
3221 // a single root file given as input
3222 // --------------------------------------
3223 }
3224 else {
3225 if (flisttxt.Contains("#"))
3226 flisttxt = flisttxt(0, flisttxt.First("#"));
3227 char *fullpath = gSystem->ConcatFileName(gSystem->DirName(flisttxt), gSystem->BaseName(flisttxt));
3228 contents->Add(new TSystemFile(fullpath, gSystem->DirName(fullpath)));// add file to the list
3229 delete fullpath;
3230 };
3231 // ---------------------------------------------------------------------------------
3232 // case 2 : no input file/file-list provided --> read all files insede the directory
3233 // ---------------------------------------------------------------------------------
3234 }
3235 else {
3236
3237 cout << "No input file list given." << endl;
3238 cout << "Check for existing root files." << endl;
3239 // cout << "Warking directory: "<< gSystem->WorkingDirectory()<< endl;
3240 if (ddir == "") {
3241 ddir = wdir;
3242 cout << "Level2 data directory : " << ddir << endl;
3243 };
3244
3245 TSystemDirectory *datadir = new TSystemDirectory(gSystem->BaseName(ddir), ddir);
3246 TList *temp = datadir->GetListOfFiles();
3247 if (!temp)
3248 return 0;
3249 // temp->Print();
3250 // cout << "*************************************" << endl;
3251
3252 TIter next(temp);
3253 TSystemFile *questo = 0;
3254
3255 while ((questo = (TSystemFile*) next())) {
3256 TString name = questo-> GetName();
3257 if (name.EndsWith(".root")) {
3258 // const char *fullpath = gSystem->FindFile(ddir,name);
3259 // char *fullpath;
3260 // gSystem->IsFileInIncludePath(name,&fullpath);
3261 char *fullpath = gSystem->ConcatFileName(gSystem->DirName(ddir), gSystem->BaseName(name));
3262 contents->Add(new TSystemFile(fullpath, gSystem->DirName(fullpath)));
3263 delete fullpath;
3264 };
3265 }
3266 delete temp;
3267 delete datadir;
3268
3269 };
3270 gSystem->ChangeDirectory(wdir); // back to the working directory
3271 // cout << endl << "Selected files:" << endl;
3272 // contents->Print();
3273 cout << contents->GetEntries() << " files \n";
3274 // cout << endl;
3275 // cout << "Working directory: "<< gSystem->WorkingDirectory()<< endl;
3276 return contents;
3277 }
3278 ;
3279 //--------------------------------------
3280 //
3281 //
3282 //--------------------------------------
3283 /**
3284 * Get the Pamela detector chains from a list of files and make them friends.
3285 * @param fl Pointer to a TList of TSystemFiles
3286 * @param detlist String to select trees to be included
3287 * @return Pointer to a TChain
3288 */
3289 TChain *PamLevel2::GetPamTree(TList *fl, TString detlist) {
3290 //
3291 //
3292 //
3293 if (pam_tree) {
3294 printf("WARNING: TChain *PamLevel2::GetPamTree(TList *fl, TString detlist) -- pam_tree already exists!\n ");
3295 return pam_tree;
3296 };
3297 //
3298
3299 TChain *Trout = 0;
3300
3301 // -------------------------------------------
3302 // set flags to include/exclude trees/branches
3303 // -------------------------------------------
3304 if (detlist.Contains("+AUTO", TString::kIgnoreCase)) {
3305 if (fl->GetEntries() > 0) {
3306 cout << "+AUTO" << endl;
3307 TFile *fprimo = new TFile(fl->At(0)->GetName());
3308 GetWhichTrees(fprimo);
3309 fprimo->Close();// AAAAARGGGGGHHHHH!!!!!!! non commentare questa riga, altrimenti si incasina il TChain
3310 fprimo->Delete();
3311 }
3312 };
3313 SetWhichTrees(detlist);
3314
3315 // -------------------------------------------
3316 // build chains
3317 // -------------------------------------------
3318 TChain *T = 0;
3319 TChain *C = 0;
3320 TChain *O = 0;
3321 TChain *R = 0;
3322 TChain *S = 0;
3323 TChain *N = 0;
3324 TChain *A = 0;
3325 TChain *B = 0;
3326 TChain *G = 0;
3327
3328 TChain *L = 0;
3329 TChain *P = 0;
3330
3331 if (TRK2 || TRK1 || TRKh)
3332 T = new TChain("Tracker");
3333 if (CAL2 || CAL1)
3334 C = new TChain("Calorimeter");
3335 if (TOF)
3336 O = new TChain("ToF");
3337 if (TRG)
3338 R = new TChain("Trigger");
3339 if (S4)
3340 S = new TChain("S4");
3341 if (ND)
3342 N = new TChain("NeutronD");
3343 if (AC)
3344 A = new TChain("Anticounter");
3345 if (ORB)
3346 B = new TChain("OrbitalInfo");
3347 if (GP)
3348 G = new TChain("h20");
3349 if (PROC)
3350 P = new TChain("ProcessingInfo");
3351 L = new TChain("SelectionList");
3352
3353 // loop over files and create chains
3354 TIter next(fl);
3355 TSystemFile *questo = 0;
3356 while ((questo = (TSystemFile*) next())) {
3357 TString name = questo->GetName();
3358 cout << "File: " << name << endl;
3359 if (CheckLevel2File(name)) {
3360 if (TRK2 || TRK1 || TRKh)
3361 T->Add(name);
3362 if (CAL1 || CAL2)
3363 C->Add(name);
3364 if (TOF)
3365 O->Add(name);
3366 if (TRG)
3367 R->Add(name);
3368 if (S4)
3369 S->Add(name);
3370 if (ND)
3371 N->Add(name);
3372 if (AC)
3373 A->Add(name);
3374 if (ORB)
3375 B->Add(name);
3376 if (GP)
3377 G->Add(name);
3378 if (P)
3379 P->Add(name);
3380 if (SELLI == 1)
3381 L->Add(name);
3382 };
3383 };
3384
3385 cout << "done data-tree chains "<< T->GetNtrees() <<" \n";
3386 cout << "----------------------------------------------------" << endl;
3387
3388 // -------------------------------------------
3389 // make friends
3390 // -------------------------------------------
3391
3392 // Tracker
3393 cout << "Friends: " << endl;
3394 if (T && (TRK2 || TRK1 || TRKh)) {
3395 if (!Trout)
3396 Trout = T;
3397 else
3398 Trout->AddFriend("Tracker");
3399 // cout << "+Tacker"<<endl;
3400 };
3401 // Calorimeter
3402 if (C && (CAL2 || CAL1)) {
3403 if (!Trout)
3404 Trout = C;
3405 else
3406 Trout->AddFriend("Calorimeter");
3407 // cout << "+Calorimeter"<<endl;
3408 };
3409 // ToF
3410 if (O && TOF) {
3411 if (!Trout)
3412 Trout = O;
3413 else
3414 Trout->AddFriend("ToF");
3415 // cout << "+ToF"<<endl;
3416 };
3417 // Trigger
3418 if (R && TRG) {
3419 if (!Trout)
3420 Trout = O;
3421 else
3422 Trout->AddFriend("Trigger");
3423 // cout << "+Trigger"<<endl;
3424 };
3425 // S4
3426 if (S && S4) {
3427 if (!Trout)
3428 Trout = O;
3429 else
3430 Trout->AddFriend("S4");
3431 // cout << "+S4"<<endl;
3432 };
3433 // Neutron Detector
3434 if (N && ND) {
3435 if (!Trout)
3436 Trout = O;
3437 else
3438 Trout->AddFriend("NeutronD");
3439 // cout << "+NeutronD"<<endl;
3440 };
3441 // Anticounters
3442 if (A && AC) {
3443 if (!Trout)
3444 Trout = O;
3445 else
3446 Trout->AddFriend("Anticounter");
3447 // cout << "+Anticounter"<<endl;
3448 };
3449 // Orbital Info
3450 if (B && ORB) {
3451 if (!Trout)
3452 Trout = O;
3453 else
3454 Trout->AddFriend("OrbitalInfo");
3455 // cout << "+OrbitalInfo"<<endl;
3456 };
3457 // GPamela
3458 if (G && GP) {
3459 if (!Trout)
3460 Trout = G;
3461 else
3462 Trout->AddFriend("h20");
3463 // cout << "+h20"<<endl;
3464 };
3465
3466 // =====================================
3467 // SET BRANCH-ADDRESS AFTER CHAIN+FRIEND
3468 // =====================================
3469 if( Trout->GetNtrees() )SetBranchAddress(Trout);
3470
3471 // ------------------------------------
3472 // finally handle selection trees...
3473 // (it is not friend of pamela tree)
3474 // ------------------------------------
3475
3476 cout << "----------------------------------------------------" << endl;
3477
3478 // Selection List
3479 if (L && SELLI == 1) {
3480 cout << ">>> Found selection-list <<<" << endl;
3481 // L->SetBranchAddress("RunEntry",&irun);
3482 L->SetBranchAddress("RunEntry", &irunt);//NEWNEW
3483 cout << "SelectionList: set branch address RunEntry" << endl;
3484 L->SetBranchAddress("EventEntry", &irunentry);
3485 cout << "SelectionList: set branch address EventEntry" << endl;
3486 /* if ( L->GetBranch("L0EventEntry") ){
3487 hasL0EE = true;
3488 L->SetBranchAddress("L0EventEntry", &il0entry);
3489 cout << "SelectionList: set branch address L0EventEntry" << endl;
3490 } else {
3491 hasL0EE = false; // backward compatibility with old preselected files...
3492 }*/
3493 sel_tree = L;
3494 // if(!Trout)Trout=O;
3495 // else Trout->AddFriend("SelectionList");
3496 cout << "+SelectionList" << endl;
3497 cout << "----------------------------------------------------" << endl;
3498 }
3499 else {
3500 // cout << "SelectionList : missing tree"<<endl;
3501 if (L)
3502 L->Delete();
3503 };
3504
3505 //ProcessingInfo EM
3506 if ( P && P->GetEntries() ){
3507 cout << "----------------------------------------------------" << endl;
3508 cout << ">>> Found ProcessingInfo <<<" << endl;
3509 // L->SetBranchAddress("RunEntry",&irun);
3510 P->SetBranchAddress("ProcInfo", &proc_obj);//NEWNEW
3511 proc_tree = P;
3512 }
3513 // --------------------------------------------
3514 // return the pamela chain with all the friends
3515 // --------------------------------------------
3516
3517 pam_tree = Trout;
3518 return Trout;
3519 }
3520
3521 //--------------------------------------
3522 //
3523 //
3524 //--------------------------------------
3525 /**
3526 * Set branch addresses for Pamela friend trees
3527 */
3528 void PamLevel2::SetBranchAddress(TTree *t) {
3529
3530 // TRK2 = TRK2 & t->GetBranchStatus("TrkLevel2");
3531 // TRK1 = TRK1 & t->GetBranchStatus("TrkLevel1");
3532 // TRKh = TRKh & t->GetBranchStatus("TrkHough");
3533 // CAL2 = CAL2 & t->GetBranchStatus("CaloLevel2");
3534 // CAL1 = CAL1 & t->GetBranchStatus("CaloLevel1");
3535 // TOF = TOF & t->GetBranchStatus("ToFLevel2");
3536 // TRG = TRG & t->GetBranchStatus("TrigLevel2");
3537 // S4 = S4 & t->GetBranchStatus("S4Level2");
3538 // ND = ND & t->GetBranchStatus("NDLevel2");
3539 // AC = AC & t->GetBranchStatus("AcLevel2");
3540 // ORB = ORB & t->GetBranchStatus("OrbitalInfo");
3541 // GP = GP & t->GetBranchStatus("h20");
3542
3543
3544 // Tracker
3545 if (TRK1) {
3546 t->SetBranchAddress("TrkLevel1", GetPointerTo("TrkLevel2"));
3547 cout << "Tracker : set branch address TrkLevel1" << endl;
3548 };
3549 if (TRK2) {
3550 t->SetBranchAddress("TrkLevel2", GetPointerTo("TrkLevel1"));
3551 cout << "Tracker : set branch address TrkLevel2" << endl;
3552 NUC = t->GetBranchStatus("TrackNuclei");
3553 EXT = t->GetBranchStatus("RecoveredTrack") && (NUC ? t->GetBranchStatus("RecoveredTrackNuclei"): true );
3554 };
3555 if (TRKh) {
3556 t->SetBranchAddress("TrkHough", GetPointerTo("TrkHough"));
3557 cout << "Tracker : set branch address TrkHough" << endl;
3558 };
3559
3560 // Calorimeter
3561 if (CAL1) {
3562 t->SetBranchAddress("CaloLevel1", GetPointerTo("CaloLevel1"));
3563 cout << "Calorimeter : set branch address CaloLevel1" << endl;
3564 };
3565 if (CAL2) {
3566 t->SetBranchAddress("CaloLevel2", GetPointerTo("CaloLevel2"));
3567 cout << "Calorimeter : set branch address CaloLevel2" << endl;
3568 };
3569
3570 // ToF
3571 if (TOF) {
3572 t->SetBranchAddress("ToFLevel2", GetPointerTo("ToFLevel2"));
3573 cout << "ToF : set branch address ToFLevel2" << endl;
3574 };
3575 // Trigger
3576 if (TRG) {
3577 t->SetBranchAddress("TrigLevel2", GetPointerTo("TrigLevel2"));
3578 cout << "Trigger : set branch address TrigLevel2" << endl;
3579 };
3580 // S4
3581 if (S4) {
3582 t->SetBranchAddress("S4Level2", GetPointerTo("S4Level2"));
3583 cout << "S4 : set branch address S4Level2" << endl;
3584 };
3585 // Neutron Detector
3586 if (ND) {
3587 t->SetBranchAddress("NDLevel2", GetPointerTo("NDLevel2"));
3588 cout << "NeutronD : set branch address NDLevel2" << endl;
3589 };
3590 // Anticounters
3591 if (AC) {
3592 t->SetBranchAddress("AcLevel2", GetPointerTo("AcLevel2"));
3593 cout << "Anticounter : set branch address AcLevel2" << endl;
3594 };
3595 // OrbitalInfo
3596 if (ORB) {
3597 t->SetBranchAddress("OrbitalInfo", GetPointerTo("OrbitalInfo"));
3598 cout << "OrbitalInfo : set branch address OrbitalInfo" << endl;
3599 };
3600 // GPamela
3601 if (GP) {
3602 // GetPointerTo("GPamela");
3603 if (!gp_obj)
3604 gp_obj = new GPamela();
3605 // gp_obj->SetBranchAddress(t); //ho dovuto fare in maniera diversa dagli altri
3606 // // t->SetBranchAddress("GPamela", GetPointerTo("GPamela"));
3607 if (SELLI)
3608 t->SetBranchAddress("GPamela", GetPointerTo("GPamela"));
3609 else
3610 gp_obj->SetBranchAddress(t); //ho dovuto fare in maniera diversa dagli altri
3611
3612 cout << "h20 : set branch address GPamela " << endl;
3613 };
3614 if(NUC){
3615
3616 cout << "Found nuclei-track branches" << endl;
3617
3618 if( !trk_nuc_obj )trk_nuc_obj = new TClonesArray("TrkTrack");
3619 if( !calo_nuc_obj)calo_nuc_obj= new TClonesArray("CaloTrkVar");
3620 if( !tof_nuc_obj)tof_nuc_obj= new TClonesArray("ToFTrkVar");
3621 if( !orb_nuc_obj)orb_nuc_obj= new TClonesArray("OrbitalInfoTrkVar");
3622 if (TRK2)t-> SetBranchAddress("TrackNuclei",&trk_nuc_obj);
3623 if (CAL2)t->GetFriend("Calorimeter")->SetBranchAddress("TrackNuclei",&calo_nuc_obj);
3624 if (TOF )t->GetFriend("ToF")-> SetBranchAddress("TrackNuclei",&tof_nuc_obj);
3625 if (ORB )t->GetFriend("OrbitalInfo")->SetBranchAddress("TrackNuclei",&orb_nuc_obj);
3626
3627 ///copy the vector content inside a "fake" object (all aother info are missing)
3628
3629 if( !trk2_nuc_obj )trk2_nuc_obj = new TrkLevel2();
3630 if( !calo2_nuc_obj )calo2_nuc_obj = new CaloLevel2();
3631 if( !tof2_nuc_obj )tof2_nuc_obj = new ToFLevel2();
3632 if( !orb2_nuc_obj )orb2_nuc_obj = new OrbitalInfo();
3633 // *(trk2_nuc_obj->GetPointerToTrackArray()) = new TClonesArray(*trk_nuc_obj);
3634 // *(calo2_nuc_obj->GetPointerToTrackArray()) = new TClonesArray(*calo_nuc_obj);
3635 // *(tof2_nuc_obj->GetPointerToTrackArray()) = new TClonesArray(*tof_nuc_obj);
3636 // *(orb2_nuc_obj->GetPointerToTrackArray()) = new TClonesArray(*orb_nuc_obj);
3637
3638 trk2_nuc_obj->SetTrackArray( trk_nuc_obj );
3639 calo2_nuc_obj->SetTrackArray( calo_nuc_obj );
3640 tof2_nuc_obj->SetTrackArray( tof_nuc_obj );
3641 orb2_nuc_obj->SetTrackArray( orb_nuc_obj );
3642
3643
3644 }
3645
3646 if(EXT){
3647
3648 cout << "Found extended tracking algorythm branches" << endl;
3649
3650 if( !trk_ext_obj )trk_ext_obj = new TClonesArray("ExtTrack");
3651 if( !calo_ext_obj)calo_ext_obj= new TClonesArray("CaloTrkVar");
3652 if( !tof_ext_obj)tof_ext_obj= new TClonesArray("ToFTrkVar");
3653 if( !orb_ext_obj)orb_ext_obj= new TClonesArray("OrbitalInfoTrkVar");
3654
3655 if (TRK2)t-> SetBranchAddress("RecoveredTrack",&trk_ext_obj);
3656 if (CAL2)t->GetFriend("Calorimeter")->SetBranchAddress("RecoveredTrack",&calo_ext_obj);
3657 if (TOF )t->GetFriend("ToF")-> SetBranchAddress("RecoveredTrack",&tof_ext_obj);
3658 if (ORB )t->GetFriend("OrbitalInfo")->SetBranchAddress("RecoveredTrack",&orb_ext_obj);
3659
3660
3661 if(NUC){
3662 if( !trk_ext_nuc_obj )trk_ext_nuc_obj = new TClonesArray("ExtTrack");
3663 if( !calo_ext_nuc_obj)calo_ext_nuc_obj= new TClonesArray("CaloTrkVar");
3664 if( !tof_ext_nuc_obj)tof_ext_nuc_obj= new TClonesArray("ToFTrkVar");
3665 if( !orb_ext_nuc_obj)orb_ext_nuc_obj= new TClonesArray("OrbitalInfoTrkVar");
3666 if (TRK2)t-> SetBranchAddress("RecoveredTrackNuclei",&trk_ext_nuc_obj);
3667 if (CAL2)t->GetFriend("Calorimeter")->SetBranchAddress("RecoveredTrackNuclei",&calo_ext_nuc_obj);
3668 if (TOF )t->GetFriend("ToF")-> SetBranchAddress("RecoveredTrackNuclei",&tof_ext_nuc_obj);
3669 if (ORB )t->GetFriend("OrbitalInfo")->SetBranchAddress("RecoveredTrackNuclei",&orb_ext_nuc_obj);
3670 }
3671 }
3672 }
3673 /**
3674 * Set branch addresses for Pamela friend trees
3675 */
3676 void PamLevel2::SetBranchAddress(TChain *t) {
3677
3678 // TRK2 = TRK2 & t->GetBranchStatus("TrkLevel2");
3679 // TRK1 = TRK1 & t->GetBranchStatus("TrkLevel1");
3680 // TRKh = TRKh & t->GetBranchStatus("TrkHough");
3681 // CAL1 = CAL1 & t->GetBranchStatus("CaloLevel1");
3682 // CAL2 = CAL2 & t->GetBranchStatus("CaloLevel2");
3683 // TOF = TOF & t->GetBranchStatus("ToFLevel2");
3684 // TRG = TRG & t->GetBranchStatus("TrigLevel2");
3685 // S4 = S4 & t->GetBranchStatus("S4Level2");
3686 // ND = ND & t->GetBranchStatus("NDLevel2");
3687 // AC = AC & t->GetBranchStatus("AcLevel2");
3688 // ORB = ORB & t->GetBranchStatus("OrbitalInfo");
3689 // GP = GP & t->GetBranchStatus("h20");
3690
3691 // Tracker
3692 if (TRK2) {
3693 t->SetBranchAddress("TrkLevel2", GetPointerTo("TrkLevel2"));
3694 cout << "Tracker : set branch address TrkLevel2" << endl;
3695 NUC = t->GetBranchStatus("TrackNuclei");
3696 EXT = t->GetBranchStatus("RecoveredTrack") && (NUC ? t->GetBranchStatus("RecoveredTrackNuclei"): true );
3697 };
3698 if (TRK1) {
3699 t->SetBranchAddress("TrkLevel1", GetPointerTo("TrkLevel1"));
3700 cout << "Tracker : set branch address TrkLevel1" << endl;
3701 };
3702 if (TRKh) {
3703 t->SetBranchAddress("TrkHough", GetPointerTo("TrkHough"));
3704 cout << "Tracker : set branch address TrkHough" << endl;
3705 };
3706
3707 // Calorimeter
3708 if (CAL2) {
3709 t->SetBranchAddress("CaloLevel2", GetPointerTo("CaloLevel2"));
3710 cout << "Calorimeter : set branch address CaloLevel2" << endl;
3711 };
3712 if (CAL1) {
3713 t->SetBranchAddress("CaloLevel1", GetPointerTo("CaloLevel1"));
3714 cout << "Calorimeter : set branch address CaloLevel1" << endl;
3715 };
3716
3717 // ToF
3718 if (TOF) {
3719 t->SetBranchAddress("ToFLevel2", GetPointerTo("ToFLevel2"));
3720 cout << "ToF : set branch address ToFLevel2" << endl;
3721 };
3722 // Trigger
3723 if (TRG) {
3724 t->SetBranchAddress("TrigLevel2", GetPointerTo("TrigLevel2"));
3725 cout << "Trigger : set branch address TrigLevel2" << endl;
3726 };
3727 // S4
3728 if (S4) {
3729 t->SetBranchAddress("S4Level2", GetPointerTo("S4Level2"));
3730 cout << "S4 : set branch address S4Level2" << endl;
3731 };
3732 // Neutron Detector
3733 if (ND) {
3734 t->SetBranchAddress("NDLevel2", GetPointerTo("NDLevel2"));
3735 cout << "NeutronD : set branch address NDLevel2" << endl;
3736 };
3737 // Anticounters
3738 if (AC) {
3739 t->SetBranchAddress("AcLevel2", GetPointerTo("AcLevel2"));
3740 cout << "Anticounter : set branch address AcLevel2" << endl;
3741 };
3742 // OrbitalInfo
3743 if (ORB) {
3744 t->SetBranchAddress("OrbitalInfo", GetPointerTo("OrbitalInfo"));
3745 cout << "OrbitalInfo : set branch address OrbitalInfo" << endl;
3746 };
3747 // GPamela
3748 // cout <<"GP "<<GP<<endl;
3749 if (GP) {
3750 // GetPointerTo("GPamela");
3751 if (!gp_obj)
3752 gp_obj = new GPamela();
3753 if (SELLI)
3754 t->SetBranchAddress("GPamela", GetPointerTo("GPamela"));
3755 else
3756 gp_obj->SetBranchAddress(t); //ho dovuto fare in maniera diversa dagli altri
3757 // gp_obj->SetBranchAddress(t); //ho dovuto fare in maniera diversa dagli altri
3758 cout << "h20 : set branch address GPamela " << endl;
3759 };
3760 // SelectionList
3761 // if(SELLI==1){
3762 // t->SetBranchAddress("RunEntry",&irunt);//NEWNEW
3763 // cout << "SelectionList: set branch address RunEntry"<<endl;
3764 // t->SetBranchAddress("EventEntry",&irunentry);
3765 // cout << "SelectionList: set branch address EventEntry"<<endl;
3766
3767 // }
3768 if(NUC){
3769
3770 cout << "Found nuclei-track branches" << endl;
3771
3772 if( !trk_nuc_obj )trk_nuc_obj = new TClonesArray("TrkTrack");
3773 if( !calo_nuc_obj)calo_nuc_obj= new TClonesArray("CaloTrkVar");
3774 if( !tof_nuc_obj)tof_nuc_obj= new TClonesArray("ToFTrkVar");
3775 if( !orb_nuc_obj)orb_nuc_obj= new TClonesArray("OrbitalInfoTrkVar");
3776 if (TRK2)t-> SetBranchAddress("TrackNuclei",&trk_nuc_obj);
3777 if (CAL2)t->GetFriend("Calorimeter")->SetBranchAddress("TrackNuclei",&calo_nuc_obj);
3778 if (TOF )t->GetFriend("ToF")-> SetBranchAddress("TrackNuclei",&tof_nuc_obj);
3779 if (ORB )t->GetFriend("OrbitalInfo")->SetBranchAddress("TrackNuclei",&orb_nuc_obj);
3780
3781 ///copy the vector content inside a "fake" object (all aother info are missing)
3782
3783 if( !trk2_nuc_obj )trk2_nuc_obj = new TrkLevel2();
3784 if( !calo2_nuc_obj )calo2_nuc_obj = new CaloLevel2();
3785 if( !tof2_nuc_obj )tof2_nuc_obj = new ToFLevel2();
3786 if( !orb2_nuc_obj )orb2_nuc_obj = new OrbitalInfo();
3787
3788 // *(trk2_nuc_obj->GetPointerToTrackArray()) = new TClonesArray(*trk_nuc_obj);
3789 // *(calo2_nuc_obj->GetPointerToTrackArray()) = new TClonesArray(*calo_nuc_obj);
3790 // *(tof2_nuc_obj->GetPointerToTrackArray()) = new TClonesArray(*tof_nuc_obj);
3791 // *(orb2_nuc_obj->GetPointerToTrackArray()) = new TClonesArray(*orb_nuc_obj);
3792 trk2_nuc_obj->SetTrackArray( trk_nuc_obj );
3793 calo2_nuc_obj->SetTrackArray( calo_nuc_obj );
3794 tof2_nuc_obj->SetTrackArray( tof_nuc_obj );
3795 orb2_nuc_obj->SetTrackArray( orb_nuc_obj );
3796
3797 }
3798 if(EXT){
3799
3800 cout << "Found extended tracking algorythm branches" << endl;
3801
3802 t->SetBranchAddress("extAlgFlag",&extAlgFlag);
3803
3804 if( !trk_ext_obj )trk_ext_obj = new TClonesArray("ExtTrack");
3805 if( !calo_ext_obj)calo_ext_obj= new TClonesArray("CaloTrkVar");
3806 if( !tof_ext_obj)tof_ext_obj= new TClonesArray("ToFTrkVar");
3807 if( !orb_ext_obj)orb_ext_obj= new TClonesArray("OrbitalInfoTrkVar");
3808
3809 if (TRK2)t-> SetBranchAddress("RecoveredTrack",&trk_ext_obj);
3810 if (CAL2)t->GetFriend("Calorimeter")->SetBranchAddress("RecoveredTrack",&calo_ext_obj);
3811 if (TOF )t->GetFriend("ToF")-> SetBranchAddress("RecoveredTrack",&tof_ext_obj);
3812 if (ORB )t->GetFriend("OrbitalInfo")->SetBranchAddress("RecoveredTrack",&orb_ext_obj);
3813
3814 if(NUC){
3815 if( !trk_ext_nuc_obj )trk_ext_nuc_obj = new TClonesArray("ExtTrack");
3816 if( !calo_ext_nuc_obj)calo_ext_nuc_obj= new TClonesArray("CaloTrkVar");
3817 if( !tof_ext_nuc_obj)tof_ext_nuc_obj= new TClonesArray("ToFTrkVar");
3818 if( !orb_ext_nuc_obj)orb_ext_nuc_obj= new TClonesArray("OrbitalInfoTrkVar");
3819 if (TRK2)t-> SetBranchAddress("RecoveredTrackNuclei",&trk_ext_nuc_obj);
3820 if (CAL2)t->GetFriend("Calorimeter")->SetBranchAddress("RecoveredTrackNuclei",&calo_ext_nuc_obj);
3821 if (TOF )t->GetFriend("ToF")-> SetBranchAddress("RecoveredTrackNuclei",&tof_ext_nuc_obj);
3822 if (ORB )t->GetFriend("OrbitalInfo")->SetBranchAddress("RecoveredTrackNuclei",&orb_ext_nuc_obj);
3823 }
3824 }
3825
3826 }
3827
3828 /**
3829 * Set the tracking algorythm
3830 * @param alg String to choose the track.
3831 * "" --> take the output of the standard tracking algorythm
3832 * "STD" --> take the output of the standard tracking algorythm
3833 * "NUC" --> take the output of the standard tracking algorythm for nuclei cluster selection
3834 * "EXT" --> in case the standard tracking algorythm has not found any track, take the output of the extended one
3835 * "EXTF" --> force the extended tracking algorythm
3836 * "NUCEXT" --> as "EXT", but for nuclei
3837 * "NUCEXTF" --> as "EXTF", but for nuclei
3838 */
3839 // void PamLevel2::SetTrackingAlgorythm(const char *alg){
3840
3841
3842 // TString s(alg);
3843 // if(s.Contains("NUC", TString::kIgnoreCase) && !NUC)
3844 // cout << "Warning! NUC algorythm requested, but branches are missing"<<endl;
3845 // if(s.Contains("EXT", TString::kIgnoreCase) && !EXT)
3846 // cout << "Warning! EXT algorythm requested, but branches are missing"<<endl;;
3847
3848 // trkAlg = alg;
3849
3850 // };
3851 // const char* PamLevel2::GetTrackingAlgorythm(){
3852
3853
3854 // cout<<endl<<" Implemented tracking algorythm: ";
3855 // cout<<endl<<" \"\" or \"STD\" --> take the output of the standard tracking algorythm";
3856 // cout<<endl<<" \"NUC\" --> take the output of the standard tracking algorythm for nuclei cluster selection";
3857 // cout<<endl<<" \"EXT\" --> in case the standard tracking algorythm has not found any track,";
3858 // cout<<endl<<" take the output of the extended one";
3859 // cout<<endl<<" \"EXTF\" --> force the extended tracking algorythm";
3860 // cout<<endl<<" \"NUCEXT\" --> as \"EXT\", but for nuclei ";
3861 // cout<<endl<<" \"NUCEXTF\" --> as \"EXTF\", but for nuclei";
3862
3863 // cout<<endl;
3864 // cout<<" <<Currently set algorythm>> "<<trkAlg<<endl;
3865 // cout<<endl;
3866
3867 // return trkAlg;
3868 // };
3869
3870
3871
3872 //--------------------------------------
3873 //
3874 //
3875 //--------------------------------------
3876 /**
3877 * Get the Run tree chain from a list of files.
3878 * @param fl Pointer to a TList of TSystemFiles
3879 * @return Pointer to a TChain
3880 */
3881 TChain *PamLevel2::GetRunTree(TList *fl) {
3882 //
3883 //
3884 //
3885 if (run_tree) {
3886 printf("WARNING: TChain *PamLevel2::GetRunTree(TList *fl) -- run_tree already exists!\n ");
3887 return run_tree;
3888 };
3889 //
3890
3891
3892 TChain *R = new TChain("Run");
3893
3894 // loop over files and create chains
3895 TIter next(fl);
3896 TSystemFile *questo = 0;
3897 while ((questo = (TSystemFile*) next())) {
3898 TString name = questo->GetName();
3899 // cout << "File: "<< name << endl;
3900 if (CheckLevel2File(name)) {
3901 R->Add(name);
3902 };
3903 }
3904
3905 cout << "done run chain "<< R->GetNtrees() <<" \n";
3906
3907
3908 if (RUN && R->GetNtrees()) {
3909
3910 R->SetBranchAddress("RunInfo", GetPointerTo("RunInfo"));
3911 cout << "Run : set branch address RunInfo" << endl;
3912 R->SetBranchAddress("SoftInfo", GetPointerTo("SoftInfo")); // Emiliano
3913 cout << "Software : set branch address SoftInfo" << endl; // Emiliano
3914
3915 irunoffset = new int[R->GetNtrees()];
3916 if (DBG) {
3917 cout << "----------------------------------------------------" << endl;
3918 cout << "irun\t | ";
3919 cout << "tree\t |";
3920 // cout << "offset\t |";
3921 cout << "RUN\t";
3922 cout << "FRAG\t";
3923 cout << "NEVENTS\t";
3924 cout << "absolute time\t\t\t";
3925 cout << "on-board time";
3926 cout << endl;
3927 }
3928 for (Int_t ii = 0; ii < R->GetEntries(); ii++) {
3929 R->GetEntry(ii);
3930 if (DBG) {
3931 cout << ii << "\t | ";
3932 cout << R->GetTreeNumber() << "\t |";
3933 // cout << R->GetChainOffset()<< "\t |";
3934 cout << GetRunInfo()->ID << "\t";
3935 cout << GetRunInfo()->ID_RUN_FRAG << "\t";
3936 cout << GetRunInfo()->NEVENTS << "\t";
3937 cout << GetRunInfo()->RUNHEADER_TIME << " <---> " << GetRunInfo()->RUNTRAILER_TIME << "\t";
3938 cout << GetRunInfo()->RUNHEADER_OBT << " <---> " << GetRunInfo()->RUNTRAILER_OBT << "\t";
3939 cout << endl;
3940 }
3941 irunoffset[R->GetTreeNumber()] = R->GetChainOffset();
3942 }
3943 cout << "N.run = " << R->GetEntries() << endl;
3944 cout << "----------------------------------------------------" << endl;
3945
3946 }
3947 // else {
3948 // delete R;
3949 // R = 0;
3950 // }
3951
3952 run_tree = R;
3953
3954 return R;
3955
3956 }
3957 //--------------------------------------
3958 //
3959 //
3960 //--------------------------------------
3961 /**
3962 * Get the Run tree from a file.
3963 * @param f Pointer to a TFile
3964 * @return Pointer to a TTree
3965 */
3966 TTree *PamLevel2::GetRunTree(TFile *f) {
3967 if (run_tree) {
3968 printf("WARNING: TTree *PamLevel2::GetRunTree(TFile *f) -- run_tree already exists!\n ");
3969 return run_tree;
3970 };
3971
3972 cout << "TTree *PamLevel2::GetRunTree(TFile *f) -- obsolte " << endl;
3973
3974 TTree *T = (TTree*) f->Get("Run");
3975
3976 if (T) {
3977 T->SetBranchAddress("RunInfo", GetPointerTo("RunInfo"));
3978 cout << "Run : set branch address RunInfo" << endl;
3979 T->SetBranchAddress("SoftInfo", GetPointerTo("SoftInfo")); // Emiliano
3980 cout << "Software : set branch address SoftInfo" << endl; // Emiliano
3981
3982 }
3983
3984 run_tree = (TChain*) T;
3985
3986 return T;
3987
3988 }
3989
3990 Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) {
3991
3992 if (DBG) printf("PamLevel2::UpdateRunInfo(Long64_t) - inside\n");
3993
3994 if (!run_tree) {
3995 cout << " Bool_t PamLevel2::UpdateRunInfo(ULong64_t iev) -- ERROR -- run tree not loaded" << endl;
3996 return false;
3997 }
3998 if (run_tree->GetEntries() <= 0) {
3999 cout << " Bool_t PamLevel2::UpdateRunInfo(ULong64_t iev) -- ERROR -- run tree is empty" << endl;
4000 return (false);
4001 }
4002
4003
4004 Int_t oldrun = irun; // store current run index
4005
4006 // -----------------------------------------------------------------------
4007 // the first time the routine is called, set run search from the beginning
4008 // -----------------------------------------------------------------------
4009
4010 if (irun < 0) {
4011 irun = 0LL;
4012 irunt = 0LL;
4013 totrunentry = 0LL;
4014 totrunentrymin = 0LL;
4015 totrunentrymax = 0LL;
4016 irunentry = 0LL;
4017 il0entry = 0LL;
4018 prevshift = 0;
4019 yprevshift = 0;
4020 prevabstime = 0;
4021 prevpktnum = 0;
4022 abstime = 0ULL;
4023 pktnum = 0;
4024 isFragment = false;
4025 run_tree->GetEntry(irun);
4026 if (!GetOrbitalInfo())
4027 cout << "** WARNING ** missing OrbitalInfo ---> run info might be not correctly updated " << endl;
4028 if ( fUseDBinRunInfo ){
4029 if (gltsync)
4030 delete gltsync; //Emiliano
4031 if (!dbc || (dbc && !dbc->IsConnected()))
4032 SetDBConnection(); //Emiliano
4033 gltsync = new GL_TIMESYNC(GetRunInfo()->ID_ROOT_L0, "ID", dbc, false); //Emiliano // the "false" means not to use level0 file (not necessary here)
4034 if (dbc){
4035 dbc->Close();// Emiliano
4036 delete dbc;
4037 dbc=0;
4038 }
4039 }
4040 }
4041 // ---------------------------------------------------------------
4042 // retrieve OBT and absolute time of the event
4043 // ---------------------------------------------------------------
4044 Long64_t obt = 0LL; // Emiliano, Long64_t GL_TIMESYNC::DBobt(UInt_t obt) since depending on the situation OBT is lowered or boosted
4045 prevabstime = abstime;
4046 prevpktnum = pktnum;
4047 if (GetOrbitalInfo()) {
4048 abstime = GetOrbitalInfo()->absTime;
4049 if ( fUseDBinRunInfo ) obt = gltsync->DBobt(GetOrbitalInfo()->OBT); // Emiliano
4050 pktnum = GetOrbitalInfo()->pkt_num; // Emiliano
4051 }
4052 else {
4053 abstime = GetRunInfo()->RUNHEADER_TIME;
4054 if ( fUseDBinRunInfo ) obt = gltsync->DBobt(GetRunInfo()->RUNHEADER_OBT); // Emiliano
4055 pktnum = GetRunInfo()->RUNHEADER_PKT; // Emiliano
4056 }
4057
4058 if (DBG){
4059 printf("0abstime %lld %lld pktnum %d %d obt %lld \n",abstime,prevabstime,pktnum,prevpktnum,obt);
4060 printf("0 rth %d %d nevents %d \n",GetRunInfo()->RUNHEADER_TIME,GetRunInfo()->RUNTRAILER_TIME,GetRunInfo()->NEVENTS);
4061 printf("0 rto %d %d \n",GetRunInfo()->RUNHEADER_OBT,GetRunInfo()->RUNTRAILER_OBT);
4062 if ( fUseDBinRunInfo ) printf("0 rto2 %lld %lld \n",gltsync->DBobt(GetRunInfo()->RUNHEADER_OBT),gltsync->DBobt(GetRunInfo()->RUNTRAILER_OBT));
4063 printf("0 bo irunentry %lld prevshift %lld irun %lld \n",irunentry,prevshift,irun);
4064 printf("0 min %lld iev %lld max %lld tot %lld \n",totrunentrymin,iev,totrunentrymax,totrunentry);
4065 }
4066
4067 // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
4068 // if it is a full file (not preselected)
4069 // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
4070 if (SELLI == 0 || SELLI == 2) { // Emiliano
4071
4072 // ---------------------------------------------------------------
4073 // increment dead and live-time counters
4074 // (only when reading a file not preselected)
4075 // ---------------------------------------------------------------
4076 if (SELLI == 0) {
4077 if (GetTrigLevel2()) {
4078 totdltime[0] += GetTrigLevel2()->dltime[0];
4079 totdltime[1] += GetTrigLevel2()->dltime[1];
4080 }
4081 totdltime[2]++;
4082 }
4083
4084 //
4085 bool a = true;
4086 bool b = true;
4087 if ( fUseDBinRunInfo ){
4088 a = false;
4089 b = false;
4090 if ( obt < gltsync->DBobt(GetRunInfo()->RUNHEADER_OBT) ) a = true;
4091 if ( obt > gltsync->DBobt(GetRunInfo()->RUNTRAILER_OBT) ) b = true;
4092 }
4093 if ( iev < totrunentrymin || iev > totrunentrymax // entry is outside run limits
4094 || iev == 0 // or it is the first entry
4095 || (!isSync && (
4096 (abstime <= GetRunInfo()->RUNHEADER_TIME && a ) // or it is outside obt limits (and abstime limits for security reasons)
4097 || (abstime >= GetRunInfo()->RUNTRAILER_TIME && b ) ))// or it is outside obt limits (and abstime limits for security reasons)
4098 ){ // check on abstime and obt needed to handle nested+DV_skipped packets
4099
4100 // check for a new run (ma prima il primo!)
4101 if (DBG){
4102 printf("1abstime %lld %lld pktnum %d %d obt %lld \n",abstime,prevabstime,pktnum,prevpktnum,obt);
4103 printf("1 rth %d %d nevents %d \n",GetRunInfo()->RUNHEADER_TIME,GetRunInfo()->RUNTRAILER_TIME,GetRunInfo()->NEVENTS);
4104 printf("1 rto %d %d \n",GetRunInfo()->RUNHEADER_OBT,GetRunInfo()->RUNTRAILER_OBT);
4105 if ( fUseDBinRunInfo ) printf("1 rto2 %lld %lld \n",gltsync->DBobt(GetRunInfo()->RUNHEADER_OBT),gltsync->DBobt(GetRunInfo()->RUNTRAILER_OBT));
4106 printf("1 bo irunentry %lld prevshift %lld irun %lld \n",irunentry,prevshift,irun);
4107 printf("1 min %lld iev %lld max %lld tot %lld \n",totrunentrymin,iev,totrunentrymax,totrunentry);
4108 }
4109 // printf("1abstime %lld %lld pktnum %d %d obt %lld \n",abstime,prevabstime,pktnum,prevpktnum,obt);
4110 // printf("1 rth %d %d nevents %d \n",GetRunInfo()->RUNHEADER_TIME,GetRunInfo()->RUNTRAILER_TIME,GetRunInfo()->NEVENTS);
4111 // printf("1 rto %d %d \n",GetRunInfo()->RUNHEADER_OBT,GetRunInfo()->RUNTRAILER_OBT);
4112 // printf("1 rto2 %lld %lld \n",gltsync->DBobt(GetRunInfo()->RUNHEADER_OBT),gltsync->DBobt(GetRunInfo()->RUNTRAILER_OBT));
4113 // printf("1 bo irunentry %lld prevshift %lld irun %lld \n",irunentry,prevshift,irun);
4114 // printf("1 min %lld iev %lld max %lld tot %lld \n",totrunentrymin,iev,totrunentrymax,totrunentry);//TOGLITOGLI
4115
4116 totrunentry = 0LL;
4117 runfirstentry = 0LL;
4118 for (Int_t r=0; r< run_tree->GetEntries();r++){
4119 // -------------------------------------------------------------------
4120 // save the index of the first entry of the run, relative to pam_tree,
4121 // and read a new run
4122 // -------------------------------------------------------------------
4123 run_tree->GetEntry(r);//update runinfo
4124 if ( r > 0 ){
4125 totrunentrymin = totrunentrymax+1;
4126 } else {
4127 totrunentrymin = 0LL;
4128 }
4129 totrunentry += GetRunInfo()->NEVENTS;
4130 totrunentrymax = totrunentry - 1 - prevshift; // prevshift is needed to handle nested+DV_skipped packets
4131 irun = r;
4132 if ( fUseDBinRunInfo ){
4133 a = false;
4134 b = false;
4135 if ( obt < gltsync->DBobt(GetRunInfo()->RUNHEADER_OBT) ) a = true;
4136 if ( obt > gltsync->DBobt(GetRunInfo()->RUNTRAILER_OBT) ) b = true;
4137 }
4138 if ( (iev >= totrunentrymin && iev <= totrunentrymax) || // entry is inside run limits
4139 ( !isSync &&
4140 ( abstime >= GetRunInfo()->RUNHEADER_TIME && a // or it is inside obt limits (and abstime limits for security reasons)
4141 && abstime <= GetRunInfo()->RUNTRAILER_TIME && b)) // or it is inside obt limits (and abstime limits for security reasons)
4142 ){ // check on abstime and obt needed to handle nested+DV_skipped packets
4143 if ( totrunentrymin > iev ){ // there is a shift (nested+DV_skipped packets)
4144 if ( !isSync ){
4145 if (DBG) printf("PamLevel2::UpdateRunInfo(Long64_t) - unconsistent iev - nevents, probable DBL0-L2 async\n");
4146 if (DBG) printf("PamLevel2::UpdateRunInfo(Long64_t) - totrunentrymin %lld iev %lld prevshift %lld totrunentrymax %lld \n",totrunentrymin,iev,prevshift,totrunentrymax);
4147 // printf("PamLevel2::UpdateRunInfo(Long64_t) - unconsistent iev - nevents, probable DBL0-L2 async\n");
4148 // printf("PamLevel2::UpdateRunInfo(Long64_t) - totrunentrymin %lld iev %lld prevshift %lld totrunentrymax %lld \n",totrunentrymin,iev,prevshift,totrunentrymax);//TOGLITOGLI
4149 prevshift += (totrunentrymin-iev); // add the new shift to total shift
4150 totrunentrymin -= (totrunentrymin-iev); // shift run position min
4151 totrunentrymax -= (totrunentrymin-iev); // shift run position max
4152 if (DBG) printf("PamLevel2::UpdateRunInfo(Long64_t) - totrunentrymin %lld iev %lld prevshift %lld totrunentrymax %lld \n",totrunentrymin,iev,prevshift,totrunentrymax);
4153 // printf("PamLevel2::UpdateRunInfo(Long64_t) - totrunentrymin %lld iev %lld prevshift %lld totrunentrymax %lld \n",totrunentrymin,iev,prevshift,totrunentrymax);//TOGLITOGLI
4154 } else {
4155 printf(" PamLevel2::UpdateRunInfo(Long64_t) ERROR! sync file but unconsistent totrunetrymin %lld and iev %lld!!! \n",totrunentrymin,iev);
4156 cout << " OK this is a bug, write to Emiliano, Emiliano.Mocchiutti@ts.infn.it " << endl;
4157 cout << "\nFor bug reporting instructions, please see for example:\n";
4158 cout << " <http://www.ts.infn.it/~mocchiut/bugs/bugs.html>.\n";
4159 cout << " " << endl;
4160 }
4161 }
4162 runfirstentry = totrunentrymin; // first entry of the run in the level2
4163
4164
4165 //
4166 if ( fUseDBinRunInfo ){
4167 if (gltsync)
4168 delete gltsync; // Emiliano
4169 if (!dbc || (dbc && !dbc->IsConnected()))
4170 SetDBConnection(); //Emiliano
4171 gltsync = new GL_TIMESYNC(GetRunInfo()->ID_ROOT_L0, "ID", dbc, false); // Emiliano
4172 TrkParams::Set(GetRunInfo(), dbc);
4173 if (dbc){
4174 dbc->Close(); // Emiliano
4175 delete dbc;
4176 dbc=0;
4177 }
4178 if (gltsync->DBobt(GetRunInfo()->RUNHEADER_OBT) > gltsync->DBobt(GetRunInfo()->RUNTRAILER_OBT)) { // Emiliano
4179 cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- WARNING -- irun " << irun
4180 << " has RUNHEADER_OBT>=RUNTRAILER_OBT " << endl;
4181 cout
4182 << " (NB!! in this case some events could be assigned to a wrong run)"
4183 << endl;
4184 }
4185 }
4186 //
4187 if (DBG) printf(" found \n");
4188 // printf(" found \n");//TOGLITOGLI
4189 //
4190 break;
4191 }
4192 } // loop over run
4193
4194 // --------------------------------------
4195 // if there was no need to update the run
4196 // ---> exit with FALSE
4197 // --------------------------------------
4198 if (irun == oldrun){
4199 if (DBG) printf(" no new run \n");
4200 // printf(" no new run \n");//TOGLITOGLI
4201 return (false);
4202 }
4203 // --------------------------------------
4204 // ... otherwise
4205 // ---> exit with TRUE
4206 // --------------------------------------
4207
4208 if (SELLI != 2) {
4209 /// decrement counters since this event belongs to a new run
4210 if (GetTrigLevel2()) {
4211 totdltime[0] -= GetTrigLevel2()->dltime[0];//live-time
4212 totdltime[1] -= GetTrigLevel2()->dltime[1];//dead-time
4213 }
4214 totdltime[2]--; //event counter
4215 if (DBG) {
4216 cout << endl;
4217 cout << "n.events : " << totdltime[2] << endl;
4218 cout << "RUN LIVE-TIME: " << totdltime[0] * 0.16 << " ms" << endl;
4219 cout << "RUN DEAD-TIME: " << totdltime[1] * 0.01 << " ms" << endl;
4220 }
4221 // add an entry
4222 if (run_tree_clone && totdltime[2] > 0)
4223 if (run_tree_clone->GetBranch("DeadLiveTime")->GetEntries() < run_tree->GetEntries())
4224 run_tree_clone->GetBranch("DeadLiveTime")->Fill();
4225 // reset counters
4226 if ( totdltime[2] > 0 ){
4227 if (GetTrigLevel2()) {
4228 totdltime[0] = GetTrigLevel2()->dltime[0];//live-time
4229 totdltime[1] = 0; //dead-time
4230 }
4231 totdltime[2] = 1; //event counter
4232 }
4233 }
4234
4235 if (DBG){
4236 cout << endl << " ))))) UPDATE RUN INFO ((((( @iev " << iev << " run " << GetRunInfo()->ID << " irun " << irun
4237 << endl;
4238 printf("2abstime %lld %lld pktnum %d %d obt %lld \n",abstime,prevabstime,pktnum,prevpktnum,obt);
4239 printf("2 rth %d %d nevents %d \n",GetRunInfo()->RUNHEADER_TIME,GetRunInfo()->RUNTRAILER_TIME,GetRunInfo()->NEVENTS);
4240 printf("2 rto %d %d \n",GetRunInfo()->RUNHEADER_OBT,GetRunInfo()->RUNTRAILER_OBT);
4241 if ( fUseDBinRunInfo ) printf("2 rto2 %lld %lld \n",gltsync->DBobt(GetRunInfo()->RUNHEADER_OBT),gltsync->DBobt(GetRunInfo()->RUNTRAILER_OBT));
4242 printf("2 bo irunentry %lld prevshift %lld irun %lld \n",irunentry,prevshift,irun);
4243 printf("2 min %lld iev %lld max %lld tot %lld \n",totrunentrymin,iev,totrunentrymax,totrunentry);
4244 }
4245 // printf("2abstime %lld %lld pktnum %d %d obt %lld \n",abstime,prevabstime,pktnum,prevpktnum,obt);
4246 // printf("2 rth %d %d nevents %d \n",GetRunInfo()->RUNHEADER_TIME,GetRunInfo()->RUNTRAILER_TIME,GetRunInfo()->NEVENTS);
4247 // printf("2 rto %d %d \n",GetRunInfo()->RUNHEADER_OBT,GetRunInfo()->RUNTRAILER_OBT);
4248 // printf("2 rto2 %lld %lld \n",gltsync->DBobt(GetRunInfo()->RUNHEADER_OBT),gltsync->DBobt(GetRunInfo()->RUNTRAILER_OBT));
4249 // printf("2 bo irunentry %lld prevshift %lld irun %lld \n",irunentry,prevshift,irun);
4250 // printf("2 min %lld iev %lld max %lld tot %lld \n",totrunentrymin,iev,totrunentrymax,totrunentry);//TOGLITOGLI
4251
4252 return (true);
4253 } // need for run upgrade
4254 if (DBG) printf("return false\n");
4255 return (false);
4256 }// SELLI = 0 SELLI = 2
4257
4258 // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
4259 // if it is a preselected file (there is SelectionList)
4260 // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
4261 // irun = run entry relative to the chain
4262 // irunt = run entry relative to the tree
4263 if (SELLI == 1) {
4264 sel_tree->GetEntry(iev);// read irunt from SelectionList
4265 irun = irunt + irunoffset[sel_tree->GetTreeNumber()];//NEWNEW
4266 if (irun != oldrun) {
4267 if (irun < run_tree->GetEntries())
4268 run_tree->GetEntry(irun);
4269 // check if the time is ok (with merged files it is not...)
4270 // if not loop over run and look for the proper entry
4271 bool SECONDO_GIRO = false;
4272 // Long64_t irun_start = irun;
4273 int offset_start = irunoffset[sel_tree->GetTreeNumber()];
4274 while (((!(abstime >= GetRunInfo()->RUNHEADER_TIME && // check on absolute time (s)
4275 abstime <= GetRunInfo()->RUNTRAILER_TIME)
4276 // ||
4277 // !(obt >= GetRunInfo()->RUNHEADER_OBT && // additional check on OBT (ms)
4278 // obt <= GetRunInfo()->RUNTRAILER_OBT)
4279 ) || GetRunInfo()->NEVENTS == 0)
4280 // && irun < run_tree->GetEntries()
4281 ) {
4282
4283 if (DBG) {
4284 cout << " (test) ";
4285 cout << " tree " << sel_tree->GetTreeNumber();
4286 cout << " irunt " << irunt;
4287 cout << " offset " << irunoffset[sel_tree->GetTreeNumber()];
4288 cout << " abs " << abstime;
4289 cout << " >> " << GetRunInfo()->RUNHEADER_TIME << " " << GetRunInfo()->RUNTRAILER_TIME;
4290 cout << " obt " << obt;
4291 cout << " >> " << GetRunInfo()->RUNHEADER_OBT << " " << GetRunInfo()->RUNTRAILER_OBT;
4292 cout << " *** JUMP RUN *** irun " << irun;
4293 cout << endl;
4294 }
4295 // irun++;
4296 irunoffset[sel_tree->GetTreeNumber()]++;
4297 irun = irunt + irunoffset[sel_tree->GetTreeNumber()];//NEWNEW
4298 if (irun == run_tree->GetEntries() && SECONDO_GIRO) {
4299 // if(irun == irun_start ){
4300 cout << " ...grrrvzzkhhhajsdkj!!!! " << endl;
4301 irunoffset[sel_tree->GetTreeNumber()] = offset_start;
4302 return false;
4303 }
4304 if (irun >= run_tree->GetEntries() || irun < 0) {
4305 cout << "irun = " << irun << " >> search from the beginning... <<" << endl;
4306 SECONDO_GIRO = true;
4307 irun = 0;
4308 irunoffset[sel_tree->GetTreeNumber()] = -irunt;
4309 }
4310 run_tree->GetEntry(irun);
4311 }
4312
4313 if (DBG) {
4314 cout << " (test) ";
4315 cout << " tree " << sel_tree->GetTreeNumber();
4316 cout << " irunt " << irunt;
4317 cout << " offset " << irunoffset[sel_tree->GetTreeNumber()];
4318 cout << " abs " << abstime;
4319 cout << " >> " << GetRunInfo()->RUNHEADER_TIME << " " << GetRunInfo()->RUNTRAILER_TIME;
4320 cout << " obt " << obt;
4321 cout << " >> " << GetRunInfo()->RUNHEADER_OBT << " " << GetRunInfo()->RUNTRAILER_OBT;
4322 }
4323 if (DBG)
4324 cout << endl;
4325 if (DBG)
4326 cout << endl << " ))))) UPDATE RUN INFO ((((( @iev " << iev << " run " << GetRunInfo()->ID << " (run-entry "
4327 << irun << ")" << endl;
4328 // ----------------------------------------------------
4329 // update the tracker parameters
4330 // (non ho trovato nessun altro modo sicuro di farlo...)
4331 // ----------------------------------------------------
4332 if ( fUseDBinRunInfo ){
4333 if (!dbc || (dbc && !dbc->IsConnected()))
4334 SetDBConnection();
4335 TrkParams::Set(GetRunInfo(), dbc);
4336 if (dbc){
4337 dbc->Close();
4338 delete dbc;
4339 dbc=0;
4340 }
4341 }
4342 // cout << endl;
4343 prevshift = 0;
4344 yprevshift = 0;
4345 return true;
4346 }
4347 return false;
4348 }
4349
4350 return false;
4351 //
4352 }
4353
4354 /**
4355 * Update the runinfo informations (to be used to have Run infos event by event basis)
4356 * @param run Pointer to the chain/tree which contains run infos
4357 * @return true if a new run has been read, false if it is still the same run
4358 */
4359 Bool_t PamLevel2::UpdateRunInfo(TChain *run, Long64_t iev) {
4360 return (UpdateRunInfo(iev));
4361 }
4362
4363 /**
4364 * Update the runinfo informations (to be used to have Run infos event by event basis)
4365 * @param run Pointer to the chain/tree which contains run infos
4366 * @return true if a new run has been read, false if it is still the same run
4367 */
4368 Bool_t PamLevel2::UpdateRunInfo(TTree *run, Long64_t iev) {
4369 return (UpdateRunInfo((TChain*) run, iev));
4370 }
4371
4372
4373 //--------------------------------------
4374 //
4375 //
4376 //--------------------------------------
4377 /**
4378 * Set which trees shoul be analysed
4379 * @param detlist TString containing the sequence of trees required
4380 */
4381 void PamLevel2::SetWhichTrees(TString detlist) {
4382
4383 // if(detlist.IsNull() || detlist.Contains("+ALL", TString::kIgnoreCase)){
4384 if (detlist.Contains("+ALL", TString::kIgnoreCase)) {
4385
4386 cout << " ======================================================== " << endl;
4387 cout << " (( WARNING )) " << endl;
4388 cout << " The meaning of the option +ALL has changed!! " << endl;
4389 cout << " Now it includes really all (level0+level1+level2+gpamela)" << endl;
4390 cout << " and the file is discarded if it does not contain " << endl;
4391 cout << " all trees or if level0 files are not available!! " << endl;
4392 cout << " ======================================================== " << endl;
4393
4394 CAL0 = true;
4395 CAL1 = true;
4396 CAL2 = true;
4397 TRK2 = true;
4398 TRK1 = true;
4399 TRKh = true;
4400 TRK0 = true;
4401 TRG = true;
4402 TOF = true;
4403 TOF0 = true;
4404 S4 = true;
4405 ND = true;
4406 AC = true;
4407 ORB = true;
4408 GP = true;
4409 }
4410 else if (detlist.Contains("-ALL", TString::kIgnoreCase)) {
4411 CAL0 = false;
4412 CAL1 = false;
4413 CAL2 = false;
4414 TRK2 = false;
4415 TRK1 = false;
4416 TRKh = false;
4417 TRK0 = false;
4418 TRG = false;
4419 TOF = false;
4420 TOF0 = false;
4421 S4 = false;
4422 ND = false;
4423 AC = false;
4424 ORB = false;
4425 GP = false;
4426 };
4427
4428 // -------------------------------------------------------------------------
4429 if (detlist.Contains("CAL1", TString::kIgnoreCase)) {
4430 if (detlist.Contains("-CAL1", TString::kIgnoreCase))
4431 CAL1 = false;
4432 if (detlist.Contains("+CAL1", TString::kIgnoreCase))
4433 CAL1 = true;
4434 };
4435
4436 if (detlist.Contains("CAL0", TString::kIgnoreCase)) {
4437 if (detlist.Contains("-CAL0", TString::kIgnoreCase))
4438 CAL0 = false;
4439 if (detlist.Contains("+CAL0", TString::kIgnoreCase))
4440 CAL0 = true;
4441 };
4442
4443 if (detlist.Contains("CAL2", TString::kIgnoreCase)) {
4444 if (detlist.Contains("-CAL2", TString::kIgnoreCase))
4445 CAL2 = false;
4446 if (detlist.Contains("+CAL2", TString::kIgnoreCase))
4447 CAL2 = true;
4448 };
4449
4450 if (detlist.Contains("+CAL", TString::kIgnoreCase) && !CAL1 && !CAL2)
4451 CAL2 = true;
4452 if (detlist.Contains("-CAL", TString::kIgnoreCase) && CAL1 && CAL2) {
4453 CAL2 = false;
4454 CAL1 = false;
4455 }
4456 // -------------------------------------------------------------------------
4457 if (detlist.Contains("TRK0", TString::kIgnoreCase)) {
4458 if (detlist.Contains("-TRK0", TString::kIgnoreCase))
4459 TRK0 = false;
4460 if (detlist.Contains("+TRK0", TString::kIgnoreCase))
4461 TRK0 = true;
4462 };
4463
4464 if (detlist.Contains("TRK1", TString::kIgnoreCase)) {
4465 if (detlist.Contains("-TRK1", TString::kIgnoreCase))
4466 TRK1 = false;
4467 if (detlist.Contains("+TRK1", TString::kIgnoreCase))
4468 TRK1 = true;
4469 };
4470
4471 if (detlist.Contains("TRK2", TString::kIgnoreCase)) {
4472 if (detlist.Contains("-TRK2", TString::kIgnoreCase))
4473 TRK2 = false;
4474 if (detlist.Contains("+TRK2", TString::kIgnoreCase))
4475 TRK2 = true;
4476 };
4477
4478 if (detlist.Contains("TRKh", TString::kIgnoreCase)) {
4479 if (detlist.Contains("-TRKh", TString::kIgnoreCase))
4480 TRKh = false;
4481 if (detlist.Contains("+TRKh", TString::kIgnoreCase))
4482 TRKh = true;
4483 };
4484
4485 if (detlist.Contains("+TRK", TString::kIgnoreCase) && !TRK1 && !TRK2 && !TRKh)
4486 TRK2 = true;
4487 if (detlist.Contains("-TRK", TString::kIgnoreCase) && TRK1 && TRK2 && TRKh) {
4488 TRK2 = false;
4489 TRK1 = false;
4490 TRKh = false;
4491 }
4492 // -------------------------------------------------------------------------
4493
4494 if (detlist.Contains("-TRG", TString::kIgnoreCase))
4495 TRG = false;
4496 else if (detlist.Contains("+TRG", TString::kIgnoreCase))
4497 TRG = true;
4498
4499 if (detlist.Contains("-TOF", TString::kIgnoreCase))
4500 TOF = false;
4501 else if (detlist.Contains("+TOF", TString::kIgnoreCase))
4502 TOF = true;
4503
4504 if (detlist.Contains("-TOF0", TString::kIgnoreCase))
4505 TOF0 = false;
4506 else if (detlist.Contains("+TOF0", TString::kIgnoreCase))
4507 TOF0 = true;
4508
4509 if (detlist.Contains("-S4", TString::kIgnoreCase))
4510 S4 = false;
4511 else if (detlist.Contains("+S4", TString::kIgnoreCase))
4512 S4 = true;
4513
4514 if (detlist.Contains("-ND", TString::kIgnoreCase))
4515 ND = false;
4516 else if (detlist.Contains("+ND", TString::kIgnoreCase))
4517 ND = true;
4518
4519 if (detlist.Contains("-AC", TString::kIgnoreCase))
4520 AC = false;
4521 else if (detlist.Contains("+AC", TString::kIgnoreCase))
4522 AC = true;
4523
4524 if (detlist.Contains("-ORB", TString::kIgnoreCase))
4525 ORB = false;
4526 else if (detlist.Contains("+ORB", TString::kIgnoreCase))
4527 ORB = true;
4528
4529 if (detlist.Contains("-GP", TString::kIgnoreCase))
4530 GP = false;
4531 else if (detlist.Contains("+GP", TString::kIgnoreCase))
4532 GP = true;
4533
4534 cout << "tree/branch list from input --> ";
4535 if (TRK0)
4536 cout << "TRK0 ";
4537 if (TRK1)
4538 cout << "TRK1 ";
4539 if (TRK2)
4540 cout << "TRK2 ";
4541 if (TRKh)
4542 cout << "TRKH ";
4543 if (CAL0)
4544 cout << "CAL0 ";
4545 if (CAL1)
4546 cout << "CAL1 ";
4547 if (CAL2)
4548 cout << "CAL2 ";
4549 if (TOF)
4550 cout << "TOF ";
4551 if (TRG)
4552 cout << "TRG ";
4553 if (AC)
4554 cout << "AC ";
4555 if (ND)
4556 cout << "ND ";
4557 if (S4)
4558 cout << "S4 ";
4559 if (ORB)
4560 cout << "ORB ";
4561 if (GP)
4562 cout << "GP ";
4563 cout << endl;
4564 // cout<< "Set detector list --> ";
4565 // if(TRK1)cout<<"TRK1 ";
4566 // if(TRK2)cout<<"TRK2 ";
4567 // if(TRKh)cout<<"TRKH ";
4568 // if(CAL1)cout<<"CAL1 ";
4569 // if(CAL2)cout<<"CAL2 ";
4570 // if(TOF0)cout<<"TOF0 ";
4571 // if(TOF)cout<<"TOF ";
4572 // if(TRG)cout<<"TRG ";
4573 // if(AC)cout<<"AC ";
4574 // if(ND)cout<<"ND ";
4575 // if(S4)cout<<"S4 ";
4576 // if(ORB)cout<<"ORB ";
4577 // cout << endl;
4578
4579 }
4580 ;
4581
4582 /**
4583 * Set tree/branch detector flags from the content of a tree
4584 */
4585 void PamLevel2::GetWhichTrees(TFile* f) {
4586
4587 // cout << "void PamLevel2::GetWhichTrees(TFile* f) --- WARNING!! --- ...potrebbe non funzionare "<<endl;
4588 // -----------
4589 // reset flags
4590 // -----------
4591 CAL1 = false;
4592 CAL2 = false;
4593 TRK2 = false;
4594 TRK1 = false;
4595 TRKh = false;
4596 TRG = false;
4597 TOF = false;
4598 S4 = false;
4599 ND = false;
4600 AC = false;
4601 ORB = false;
4602 GP = false;
4603
4604 RUN = false;
4605
4606 // cout << "Checking file: "<<f->GetName()<<endl;
4607 if (!f || f->IsZombie()) {
4608 cout << "File: " << f->GetName() << " Non valid root file" << endl;
4609 fDiscarded = true;
4610 return;
4611 }
4612
4613 TList *lk = f->GetListOfKeys();
4614 if (!lk)
4615 return;
4616 TIter next(lk);
4617 TKey *key = 0;
4618
4619 Int_t nev = 0;
4620
4621 while ((key = (TKey*) next())) {
4622
4623 if (!strcmp(key->GetName(), "Run"))
4624 RUN = true;
4625
4626 //=========================================================
4627 if (!strcmp(key->GetName(), "Trigger")) {
4628 TRG = true;
4629 Int_t nevt = ((TTree*) f->Get("Trigger"))->GetEntries();
4630 if (nev && nevt != nev) {
4631 cout << "File: " << f->GetName() << " Trigger tree has " << nevt << " events instead of " << nev << endl;
4632 TRG = false;
4633 }
4634 else
4635 nev = nevt;
4636 }
4637 //=========================================================
4638 if (!strcmp(key->GetName(), "ToF")) {
4639 TOF = true;
4640 Int_t nevt = ((TTree*) f->Get("ToF"))->GetEntries();
4641 if (nev && nevt != nev) {
4642 cout << "File: " << f->GetName() << " ToF tree has " << nevt << " events instead of " << nev << endl;
4643 TOF = false;
4644 }
4645 else
4646 nev = nevt;
4647 }
4648 //=========================================================
4649 if (!strcmp(key->GetName(), "S4")) {
4650 S4 = true;
4651 Int_t nevt = ((TTree*) f->Get("S4"))->GetEntries();
4652 if (nev && nevt != nev) {
4653 cout << "File: " << f->GetName() << " S4 tree has " << nevt << " events instead of " << nev << endl;
4654 S4 = false;
4655 }
4656 else
4657 nev = nevt;
4658 }
4659 //=========================================================
4660
4661 if (!strcmp(key->GetName(), "NeutronD")) {
4662 ND = true;
4663 Int_t nevt = ((TTree*) f->Get("NeutronD"))->GetEntries();
4664 if (nev && nevt != nev) {
4665 cout << "File: " << f->GetName() << "NeutronD tree has " << nevt << " events instead of " << nev << endl;
4666 ND = false;
4667 }
4668 else
4669 nev = nevt;
4670 }
4671 //=========================================================
4672 if (!strcmp(key->GetName(), "Anticounter")) {
4673 AC = true;
4674 Int_t nevt = ((TTree*) f->Get("Anticounter"))->GetEntries();
4675 if (nev && nevt != nev) {
4676 cout << "File: " << f->GetName() << " Anticounter tree has " << nevt << " events instead of " << nev << endl;
4677 AC = false;
4678 }
4679 else
4680 nev = nevt;
4681 }
4682 //=========================================================
4683 if (!strcmp(key->GetName(), "OrbitalInfo")) {
4684 ORB = true;
4685 Int_t nevt = ((TTree*) f->Get("OrbitalInfo"))->GetEntries();
4686 if (nev && nevt != nev) {
4687 cout << "File: " << f->GetName() << " OrbitalInfo tree has " << nevt << " events instead of " << nev << endl;
4688 ORB = false;
4689 }
4690 else
4691 nev = nevt;
4692 }
4693 //=========================================================
4694 if (!strcmp(key->GetName(), "Tracker")) {
4695 TTree *T = (TTree*) f->Get("Tracker");
4696 for (Int_t i = 0; i < T->GetListOfBranches()->GetEntries(); i++) {
4697 TString name = T->GetListOfBranches()->At(i)->GetName();
4698 if (!name.CompareTo("TrkLevel1"))
4699 TRK1 = true;
4700 if (!name.CompareTo("TrkLevel2"))
4701 TRK2 = true;
4702 if (!name.CompareTo("TrkHough"))
4703 TRKh = true;
4704 };
4705 Int_t nevt = T->GetEntries();
4706 if (nev && nevt != nev) {
4707 cout << "File: " << f->GetName() << " Tracker tree has " << nevt << " events instead of " << nev << endl;
4708 TRK1 = false;
4709 TRK2 = false;
4710 TRKh = false;
4711 }
4712 else
4713 nev = nevt;
4714 // T->Delete();
4715 };
4716 //=========================================================
4717 if (!strcmp(key->GetName(), "Calorimeter")) {
4718 TTree *T = (TTree*) f->Get("Calorimeter");
4719 for (Int_t i = 0; i < T->GetListOfBranches()->GetEntries(); i++) {
4720 TString name = T->GetListOfBranches()->At(i)->GetName();
4721 if (!name.CompareTo("CaloLevel1"))
4722 CAL1 = true;
4723 if (!name.CompareTo("CaloLevel2"))
4724 CAL2 = true;
4725 };
4726 Int_t nevt = T->GetEntries();
4727 if (nev && nevt != nev) {
4728 cout << "File: " << f->GetName() << " Calorimeter tree has " << nevt << " events instead of " << nev << endl;
4729 CAL1 = false;
4730 CAL2 = false;
4731 }
4732 else
4733 nev = nevt;
4734 // T->Delete();
4735 };
4736 //=========================================================
4737 if (!strcmp(key->GetName(), "h20")) {
4738 GP = true;
4739 Int_t nevt = ((TTree*) f->Get("h20"))->GetEntries();
4740 if (nev && nevt != nev) {
4741 cout << "File: " << f->GetName() << " h20 tree has " << nevt << " events instead of " << nev << endl;
4742 GP = false;
4743 }
4744 else
4745 nev = nevt;
4746 }
4747
4748 };
4749
4750 // delete lk;
4751
4752 cout << "tree/branch list from file --> ";
4753 if (TRK1)
4754 cout << "TRK1 ";
4755 if (TRK2)
4756 cout << "TRK2 ";
4757 if (TRKh)
4758 cout << "TRKH ";
4759 if (CAL1)
4760 cout << "CAL1 ";
4761 if (CAL2)
4762 cout << "CAL2 ";
4763 if (TOF)
4764 cout << "TOF ";
4765 if (TRG)
4766 cout << "TRG ";
4767 if (AC)
4768 cout << "AC ";
4769 if (ND)
4770 cout << "ND ";
4771 if (S4)
4772 cout << "S4 ";
4773 if (ORB)
4774 cout << "ORB ";
4775 if (GP)
4776 cout << "GP ";
4777 cout << endl;
4778
4779 return;
4780
4781 }
4782 ;
4783
4784 //--------------------------------------
4785 //
4786 //
4787 //--------------------------------------
4788 /**
4789 * Check if a file contains selected Pamela Level2 trees.
4790 * @param name File name
4791 * @return true if the file is ok.
4792 */
4793 Bool_t PamLevel2::CheckLevel2File(TString name) {
4794
4795 Bool_t CAL1__ok = false;
4796 Bool_t CAL2__ok = false;
4797 Bool_t TRK2__ok = false;
4798 Bool_t TRK1__ok = false;
4799 Bool_t TRKh__ok = false;
4800 Bool_t TRG__ok = false;
4801 Bool_t TOF__ok = false;
4802 Bool_t S4__ok = false;
4803 Bool_t ND__ok = false;
4804 Bool_t AC__ok = false;
4805 Bool_t ORB__ok = false;
4806 Bool_t GP__ok = false;
4807
4808 Bool_t RUN__ok = false;
4809
4810 Bool_t SELLI__ok = false;
4811
4812 // cout << "Checking file: "<<name<<endl;
4813 TFile *f = new TFile(name.Data());
4814 if (!f || f->IsZombie()) {
4815 cout << "File: " << f->GetName() << " discarded ---- Non valid root file" << endl;
4816 fDiscarded = true;
4817 return false;
4818 }
4819 // cout << "Get list of keys: "<<f<<endl;
4820 TList *lk = f->GetListOfKeys();
4821 // lk->Print();
4822 TIter next(lk);
4823 TKey *key = 0;
4824
4825 Int_t nev = 0;
4826
4827 while ((key = (TKey*) next())) {
4828
4829 // cout << key->GetName() << endl;
4830 // cout << key->GetName() << ""<<key->GetClassName()<<endl;
4831 // cout << " Get tree: " << f->Get(key->GetName())<<endl;
4832 // nev_previous = nev;
4833 // cout << " n.entries "<< nev <<endl;
4834 // if( key->GetClassName()=="TTree" && nev_previous && nev != nev_previous ){
4835 // nev = ((TTree*)f->Get(key->GetName()))->GetEntries();
4836 // cout << "File: "<< f->GetName() <<" discarded ---- "<< key->GetName() << " tree: n.entries does not match "<<nev<<" "<<nev_previous<< endl;
4837 // return false;
4838 // };
4839
4840 //=========================================================
4841 // check if the file
4842
4843
4844 if (!strcmp(key->GetName(), "Run"))
4845 RUN__ok = true;
4846
4847 //=========================================================
4848 if (!strcmp(key->GetName(), "SelectionList")) {
4849 SELLI__ok = true;
4850 if (SELLI == 1) {
4851 Int_t nevt = ((TTree*) f->Get("SelectionList"))->GetEntries();
4852 if (nev && nevt != nev) {
4853 cout << "File: " << f->GetName() << " discarded ---- SelectionList tree has " << nevt
4854 << " events instead of " << nev << endl;
4855 fDiscarded = true;
4856 return false;
4857 }
4858 nev = nevt;
4859 }
4860 }
4861
4862 //=========================================================
4863 if (!strcmp(key->GetName(), "Trigger")) {
4864 TRG__ok = true;
4865 if (TRG) {
4866 Int_t nevt = ((TTree*) f->Get("Trigger"))->GetEntries();
4867 if (nev && nevt != nev) {
4868 cout << "File: " << f->GetName() << " discarded ---- Trigger tree has " << nevt << " events instead of "
4869 << nev << endl;
4870 fDiscarded = true;
4871 return false;
4872 }
4873 nev = nevt;
4874 }
4875 }
4876 //=========================================================
4877 if (!strcmp(key->GetName(), "ToF")) {
4878 TOF__ok = true;
4879 if (TOF) {
4880 Int_t nevt = ((TTree*) f->Get("ToF"))->GetEntries();
4881 if (nev && nevt != nev) {
4882 cout << "File: " << f->GetName() << " discarded ---- ToF tree has " << nevt << " events instead of " << nev
4883 << endl;
4884 fDiscarded = true;
4885 return false;
4886 }
4887 nev = nevt;
4888 }
4889 }
4890 //=========================================================
4891 if (!strcmp(key->GetName(), "S4")) {
4892 S4__ok = true;
4893 if (S4) {
4894 Int_t nevt = ((TTree*) f->Get("S4"))->GetEntries();
4895 if (nev && nevt != nev) {
4896 cout << "File: " << f->GetName() << " discarded ---- S4 tree has " << nevt << " events instead of " << nev
4897 << endl;
4898 fDiscarded = true;
4899 return false;
4900 }
4901 nev = nevt;
4902 }
4903 }
4904 //=========================================================
4905
4906 if (!strcmp(key->GetName(), "NeutronD")) {
4907 ND__ok = true;
4908 if (ND) {
4909 Int_t nevt = ((TTree*) f->Get("NeutronD"))->GetEntries();
4910 if (nev && nevt != nev) {
4911 cout << "File: " << f->GetName() << " discarded ---- NeutronD tree has " << nevt << " events instead of "
4912 << nev << endl;
4913 fDiscarded = true;
4914 return false;
4915 }
4916 nev = nevt;
4917 }
4918 }
4919 //=========================================================
4920 if (!strcmp(key->GetName(), "Anticounter")) {
4921 AC__ok = true;
4922 if (AC) {
4923 Int_t nevt = ((TTree*) f->Get("Anticounter"))->GetEntries();
4924 if (nev && nevt != nev) {
4925 cout << "File: " << f->GetName() << " discarded ---- Anticounter tree has " << nevt << " events instead of "
4926 << nev << endl;
4927 fDiscarded = true;
4928 return false;
4929 }
4930 nev = nevt;
4931 }
4932 }
4933 //=========================================================
4934 if (!strcmp(key->GetName(), "OrbitalInfo")) {
4935 ORB__ok = true;
4936 if (ORB) {
4937 Int_t nevt = ((TTree*) f->Get("OrbitalInfo"))->GetEntries();
4938 if (nev && nevt != nev) {
4939 cout << "File: " << f->GetName() << " discarded ---- OrbitalInfo tree has " << nevt << " events instead of "
4940 << nev << endl;
4941 fDiscarded = true;
4942 return false;
4943 }
4944 nev = nevt;
4945 }
4946 }
4947 //=========================================================
4948 if (!strcmp(key->GetName(), "Tracker")) {
4949 TTree *T = (TTree*) f->Get("Tracker");
4950 if (TRK1 || TRK2 || TRKh) {
4951 Int_t nevt = T->GetEntries();
4952 if (nev && nevt != nev) {
4953 cout << "File: " << f->GetName() << " discarded ---- Tracker tree has " << nevt << " events instead of "
4954 << nev << endl;
4955 fDiscarded = true;
4956 return false;
4957 }
4958 nev = nevt;
4959 }
4960 for (Int_t i = 0; i < T->GetListOfBranches()->GetEntries(); i++) {
4961 TString name = T->GetListOfBranches()->At(i)->GetName();
4962 if (!name.CompareTo("TrkLevel1"))
4963 TRK1__ok = true;
4964 if (!name.CompareTo("TrkLevel2"))
4965 TRK2__ok = true;
4966 if (!name.CompareTo("TrkHough"))
4967 TRKh__ok = true;
4968 };
4969 T->Delete();
4970 };
4971 //=========================================================
4972 if (!strcmp(key->GetName(), "Calorimeter")) {
4973 TTree *T = (TTree*) f->Get("Calorimeter");
4974 if (CAL1 || CAL2) {
4975 Int_t nevt = T->GetEntries();
4976 if (nev && nevt != nev) {
4977 cout << "File: " << f->GetName() << " discarded ---- Calorimeter tree has " << nevt << " events instead of "
4978 << nev << endl;
4979 fDiscarded = true;
4980 return false;
4981 }
4982 nev = nevt;
4983 }
4984 for (Int_t i = 0; i < T->GetListOfBranches()->GetEntries(); i++) {
4985 TString name = T->GetListOfBranches()->At(i)->GetName();
4986 if (!name.CompareTo("CaloLevel1"))
4987 CAL1__ok = true;
4988 if (!name.CompareTo("CaloLevel2"))
4989 CAL2__ok = true;
4990 };
4991 T->Delete();
4992 };
4993 //=========================================================
4994 if (!strcmp(key->GetName(), "h20")) {
4995 ISGP = true;
4996 GP__ok = true;
4997 if (GP) {
4998 Int_t nevt = ((TTree*) f->Get("h20"))->GetEntries();
4999 if (nev && nevt != nev) {
5000 cout << "File: " << f->GetName() << " discarded ---- h20 tree has " << nevt << " events instead of " << nev
5001 << endl;
5002 fDiscarded = true;
5003 return false;
5004 }
5005 nev = nevt;
5006 }
5007 }
5008
5009 };
5010
5011 if (SELLI == -1)
5012 SELLI = (Int_t) SELLI__ok;
5013 if (SELLI == 0 && SELLI__ok) {
5014 cout << "File: " << f->GetName() << " discarded ---- found SelectionList (it is not a full-event file)" << endl;
5015 fDiscarded = true;
5016 return false;
5017 }
5018 if (SELLI == 1 && !SELLI__ok) {
5019 cout << "File: " << f->GetName() << " discarded ---- SelectionList missing" << endl;
5020 fDiscarded = true;
5021 return false;
5022 }
5023
5024 // cout << "SELLI "<<SELLI<<endl;
5025
5026 // cout<< "CheckLevel2File(TString): detector list --> ";
5027 // if(TRK1__ok)cout<<"TRK1 ";
5028 // if(TRK2__ok)cout<<"TRK2 ";
5029 // if(TRKh__ok)cout<<"TRKH ";
5030 // if(CAL1__ok)cout<<"CAL1 ";
5031 // if(CAL2__ok)cout<<"CAL2 ";
5032 // if(TOF__ok)cout<<"TOF ";
5033 // if(TRG__ok)cout<<"TRG ";
5034 // if(AC__ok)cout<<"AC ";
5035 // if(ND__ok)cout<<"ND ";
5036 // if(S4__ok)cout<<"S4 ";
5037 // if(ORB__ok)cout<<"ORB ";
5038 // cout << endl;
5039
5040
5041 if (TRK2 && TRK1__ok)
5042 TRK1 = 1;
5043 // ----------------------------------------------------------------------------
5044 // NOTA
5045 // se c'e` il level1, lo devo necessarimente leggere.
5046 // infatti (non ho capito perche`) i cluster vengono letti e allocati in memoria
5047 // comunque, ma non vengono disallocati da PamLevel2::Clear()
5048 // ----------------------------------------------------------------------------
5049
5050
5051 if (!RUN__ok) {
5052 cout << "File: " << f->GetName() << " *WARNING* ---- Missing RunInfo tree (NB: RUN infos will not be updated)"
5053 << endl;
5054 RUN = false;
5055 };
5056
5057 if (CAL1 && !CAL1__ok) {
5058 cout << "File: " << f->GetName() << " discarded ---- Missing CaloLevel1 branch" << endl;
5059 fDiscarded = true;
5060 return false;
5061 };
5062 if (CAL2 && !CAL2__ok) {
5063 cout << "File: " << f->GetName() << " discarded ---- Missing CaloLevel2 branch" << endl;
5064 fDiscarded = true;
5065 return false;
5066 };
5067 if (TRK2 && !TRK2__ok) {
5068 cout << "File: " << f->GetName() << " discarded ---- Missing TrkLevel2 branch" << endl;
5069 fDiscarded = true;
5070 return false;
5071 };
5072 if (TRK1 && !TRK1__ok) {
5073 cout << "File: " << f->GetName() << " discarded ---- Missing TrkLevel1 branch" << endl;
5074 fDiscarded = true;
5075 return false;
5076 };
5077 if (TRKh && !TRKh__ok) {
5078 cout << "File: " << f->GetName() << " discarded ---- Missing TrkHough branch" << endl;
5079 fDiscarded = true;
5080 return false;
5081 };
5082 if (ORB && !ORB__ok) {
5083 cout << "File: " << f->GetName() << " discarded ---- Missing ORB tree" << endl;
5084 fDiscarded = true;
5085 return false;
5086 };
5087 if (AC && !AC__ok) {
5088 cout << "File: " << f->GetName() << " discarded ---- Missing AC tree" << endl;
5089 fDiscarded = true;
5090 return false;
5091 };
5092 if (S4 && !S4__ok) {
5093 cout << "File: " << f->GetName() << " discarded ---- Missing S4 tree" << endl;
5094 fDiscarded = true;
5095 return false;
5096 };
5097 if (TOF && !TOF__ok) {
5098 cout << "File: " << f->GetName() << " discarded ---- Missing ToF tree" << endl;
5099 fDiscarded = true;
5100 return false;
5101 };
5102
5103 if (ND && !ND__ok) {
5104 cout << "File: " << f->GetName() << " discarded ---- Missing ND tree" << endl;
5105 fDiscarded = true;
5106 return false;
5107 };
5108 if (TRG && !TRG__ok) {
5109 cout << "File: " << f->GetName() << " discarded ---- Missing Trigger tree" << endl;
5110 fDiscarded = true;
5111 return false;
5112 };
5113 if (GP && !GP__ok) {
5114 cout << "File: " << f->GetName() << " discarded ---- Missing h20 tree" << endl;
5115 fDiscarded = true;
5116 return false;
5117 };
5118
5119 // lk->Delete();
5120 // delete lk;
5121 f->Close();
5122
5123 // cout<< "CheckLevel2File(TString): detector list --> ";
5124 // if(TRK1)cout<<"TRK1 ";
5125 // if(TRK2)cout<<"TRK2 ";
5126 // if(TRKh)cout<<"TRKH ";
5127 // if(CAL1)cout<<"CAL1 ";
5128 // if(CAL2)cout<<"CAL2 ";
5129 // if(TOF)cout<<"TOF ";
5130 // if(TRG)cout<<"TRG ";
5131 // if(AC)cout<<"AC ";
5132 // if(ND)cout<<"ND ";
5133 // if(S4)cout<<"S4 ";
5134 // if(ORB)cout<<"ORB ";
5135 // if(GP)cout<<"GP ";
5136 // cout << endl;
5137
5138 return true;
5139
5140 }
5141 ;
5142
5143 /**
5144 * Create clone-trees
5145 */
5146 void PamLevel2::CreateCloneTrees0(TChain *fChain, TFile *ofile) {
5147
5148 cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" << endl;
5149 cout << "Create clones of PAMELA trees " << endl;
5150
5151 Int_t i = 0;
5152 pam_tree_clone[i] = fChain->GetTree()->CloneTree(0);
5153 TString name = pam_tree_clone[i]->GetName();
5154 name.Append("_clone");
5155 // pam_tree_clone[i]->SetName(name.Data());
5156 cout << pam_tree_clone[i]->GetName() << endl;
5157 i++;
5158
5159 TList *li = fChain->GetListOfFriends();
5160 TIter next(li);
5161 TFriendElement* T_friend = 0;
5162 ofile->cd();
5163 while ((T_friend = (TFriendElement*) next())) {
5164 // cout<<T_friend->IsA()->GetName()<<" "<<T_friend->GetName()<<hex << T_friend->GetTree() << dec<<endl;
5165 // cout<<T_friend->GetTree()->GetName()<< endl;
5166 pam_tree_clone[i] = T_friend->GetTree()->CloneTree(0);
5167 pam_tree_clone[i]->SetAutoSave(1000000);
5168 name = pam_tree_clone[i]->GetName();
5169 name.Append("_clone");
5170 // pam_tree_clone[i]->SetName(name.Data());
5171 cout << pam_tree_clone[i]->GetName() << endl;
5172 i++;
5173 }
5174
5175 delete li;
5176
5177 cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" << endl;
5178
5179 }
5180
5181 /**
5182 * Create clone-trees
5183 */
5184 void PamLevel2::CreateCloneTrees(TFile *ofile) {
5185
5186 // if the pointer is null, create a default file
5187 if (!run_tree)
5188 return;
5189
5190 if (!ofile) {
5191 cout << "void PamLevel2::CreateCloneTrees(TFile*) -- WARNING -- Creating file: clone-tree.root " << endl;
5192 ofile = new TFile("clone-tree.root", "recreate");
5193 }
5194
5195 ofile->cd();
5196
5197 cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" << endl;
5198 cout << "Create new PAMELA trees " << endl;
5199
5200 run_tree_clone = new TTree("Run", "PAMELA Level2 data from the GL_RUN table ");
5201 run_tree_clone->Branch("RunInfo", "GL_RUN", GetPointerTo("RunInfo"));
5202 cout << "Run : branch RunInfo" << endl;
5203 run_tree_clone->Branch("SoftInfo", "SoftInfo", GetPointerTo("SoftInfo"));
5204 cout << "Run : branch SoftInfo" << endl;
5205 // ------------------
5206 // replicate run tree
5207 // ------------------
5208 // cout << "----------------------------------------------------"<<endl;
5209 // cout << "irun\t | RUN\t NEVENTS\t absolute time"<<endl;
5210 for (Int_t i = 0; i < run_tree->GetEntries(); i++) {
5211 run_tree->GetEntry(i);
5212 // cout << i<< "\t | "<<GetRunInfo()->ID<<"\t "<<GetRunInfo()->NEVENTS<< "\t "<<GetRunInfo()->RUNHEADER_TIME<<" <---> "<<GetRunInfo()->RUNTRAILER_TIME<<endl;
5213 run_tree_clone->Fill();
5214 }
5215 //
5216 // replicate processinginfo tree
5217 //
5218 if ( PROC && false ){ // EMEMEM
5219 proc_tree_clone = new TTree("ProcessingInfo","Log of data processing");
5220 proc_tree_clone->Branch("ProcInfo", "ProcInfo", GetPointerTo("ProcInfo"));
5221 cout << "ProcessingInfo: branch ProcessingInfo" << endl;
5222 // ------------------
5223 // replicate processinginfo tree
5224 // ------------------
5225 // cout << "----------------------------------------------------"<<endl;
5226 // cout << "irun\t | RUN\t NEVENTS\t absolute time"<<endl;
5227 for (Int_t i = 0; i < proc_tree->GetEntries(); i++) {
5228 proc_tree->GetEntry(i);
5229 cout << i<< "\t | "<<endl;
5230 proc_tree_clone->Fill();
5231 }
5232 if ( SELLI != 2 || true ){
5233 cout << "|| "<<endl;
5234 proc_obj->runID = 0;
5235 cout << "|||d "<<endl;
5236 TTimeStamp *dt = new TTimeStamp();
5237 proc_obj->date = dt->AsString();
5238 delete dt;
5239 cout << "|||f "<<endl;
5240 proc_obj->commandLine = Form("PamelaLevel2 was called");
5241 cout << "|||g "<<endl;
5242 proc_obj->outputFilename = "";
5243 cout << "|||h "<<endl;
5244 proc_obj->localDir = gSystem->WorkingDirectory();
5245 cout << "|||j "<<endl;
5246 proc_obj->uname = gSystem->GetFromPipe("uname -a");
5247 cout << "|||s "<<endl;
5248 proc_obj->DB = Form("mysql://%s/%s",dbc->GetHost(),dbc->GetDB());
5249 cout << "||| "<<endl;
5250 proc_tree_clone->Fill();
5251 }
5252 cout << "----------------------------------------------------" << endl;
5253 }
5254 // ------------------------------------
5255 // add branch with dead and live times
5256 // ------------------------------------
5257 if (SELLI != 2) { // EMILIANO
5258 run_tree_clone->Branch("DeadLiveTime", totdltime, "dltime[3]/l");
5259 cout << "Run : branch DeadLiveTime" << endl;
5260
5261 sel_tree_clone = new TTree("SelectionList", "List of selected events ");
5262 // sel_tree_clone->Branch("RunEntry",&irun,"runentry/L");
5263 sel_tree_clone->Branch("RunEntry", &irunt, "runentry/L");//NEWNEW
5264 sel_tree_clone->Branch("EventEntry", &irunentry, "eventry/L");
5265 // if ( hasL0EE ) sel_tree_clone->Branch("L0EventEntry", &il0entry, "l0eventry/L");
5266 };
5267
5268 Int_t i = 0;
5269 if (TRK1 || TRK2 || TRKh) {
5270 pam_tree_clone[i] = new TTree("Tracker", "PAMELA tracker level2 data ");
5271 if (TRK1) {
5272 pam_tree_clone[i]->Branch("TrkLevel1", "TrkLevel1", GetPointerTo("TrkLevel1"));
5273 pam_tree_clone[i]->BranchRef();
5274 cout << "Tracker : branch TrkLevel1" << endl;
5275 // cout << "CreateCloneTrees " << GetTrkLevel1()<<endl;
5276 };
5277 if (TRK2) {
5278 pam_tree_clone[i]->Branch("TrkLevel2", "TrkLevel2", GetPointerTo("TrkLevel2"));
5279 cout << "Tracker : branch TrkLevel2" << endl;
5280 };
5281 if (TRKh) {
5282 pam_tree_clone[i]->Branch("TrkHough", "TrkHough", GetPointerTo("TrkHough"));
5283 cout << "Tracker : branch TrkHough" << endl;
5284 };
5285 if(NUC){
5286 pam_tree_clone[i]->Branch("TrackNuclei","TClonesArray",&trk_nuc_obj);
5287 cout << "Tracker : branch TrackNuclei" << endl;
5288 }
5289 if(EXT){
5290 pam_tree_clone[i]->Branch("extAlgFlag",&extAlgFlag,"extAlgFlag/I");
5291 pam_tree_clone[i]->Branch("RecoveredTrack","TClonesArray",&trk_ext_obj);
5292 cout << "Tracker : branch RecoveredTrack" << endl;
5293 if(NUC){
5294 pam_tree_clone[i]->Branch("RecoveredTrackNuclei","TClonesArray",&trk_ext_nuc_obj);
5295 cout << "Tracker : branch RecoveredTrackNuclei" << endl;
5296 }
5297 }
5298
5299 i++;
5300 }
5301
5302 // Calorimeter
5303 if (CAL1 || CAL2) {
5304 pam_tree_clone[i] = new TTree("Calorimeter", "PAMELA calorimeter level2 data ");
5305 if (CAL1) {
5306 pam_tree_clone[i]->Branch("CaloLevel1", "CaloLevel1", GetPointerTo("CaloLevel1"));
5307 cout << "Calorimeter : branch CaloLevel1" << endl;
5308 };
5309 if (CAL2) {
5310 pam_tree_clone[i]->Branch("CaloLevel2", "CaloLevel2", GetPointerTo("CaloLevel2"));
5311 cout << "Calorimeter : branch CaloLevel2" << endl;
5312 };
5313 if(NUC){
5314 pam_tree_clone[i]->Branch("TrackNuclei","TClonesArray",&calo_nuc_obj);
5315 cout << "Calorimeter : branch TrackNuclei" << endl;
5316 }
5317 if(EXT){
5318 pam_tree_clone[i]->Branch("RecoveredTrack","TClonesArray",&calo_ext_obj);
5319 cout << "Calorimeter : branch RecoveredTrack" << endl;
5320 if(NUC){
5321 pam_tree_clone[i]->Branch("RecoveredTrackNuclei","TClonesArray",&calo_ext_nuc_obj);
5322 cout << "Calorimeter : branch RecoveredTrackNuclei" << endl;
5323 }
5324 }
5325 i++;
5326 }
5327
5328 // ToF
5329 if (TOF) {
5330 pam_tree_clone[i] = new TTree("ToF", "PAMELA ToF level2 data ");
5331 pam_tree_clone[i]->Branch("ToFLevel2", "ToFLevel2", GetPointerTo("ToFLevel2"));
5332 cout << "ToF : branch ToFLevel2" << endl;
5333 if(NUC){
5334 pam_tree_clone[i]->Branch("TrackNuclei","TClonesArray",&tof_nuc_obj);
5335 cout << "ToF : branch TrackNuclei" << endl;
5336 }
5337 if(EXT){
5338 pam_tree_clone[i]->Branch("RecoveredTrack","TClonesArray",&tof_ext_obj);
5339 cout << "ToF : branch RecoveredTrack" << endl;
5340 if(NUC){
5341 pam_tree_clone[i]->Branch("RecoveredTrackNuclei","TClonesArray",&tof_ext_nuc_obj);
5342 cout << "ToF : branch RecoveredTrackNuclei" << endl;
5343 }
5344 }
5345 i++;
5346 };
5347 // Trigger
5348 if (TRG) {
5349 pam_tree_clone[i] = new TTree("Trigger", "PAMELA trigger level2 data ");
5350 pam_tree_clone[i]->Branch("TrigLevel2", "TrigLevel2", GetPointerTo("TrigLevel2"));
5351 cout << "Trigger : branch TrigLevel2" << endl;
5352 i++;
5353 };
5354 // S4
5355 if (S4) {
5356 pam_tree_clone[i] = new TTree("S4", "PAMELA S4 level2 data ");
5357 pam_tree_clone[i]->Branch("S4Level2", "S4Level2", GetPointerTo("S4Level2"));
5358 cout << "S4 : branch S4Level2" << endl;
5359 i++;
5360 };
5361 // Neutron Detector
5362 if (ND) {
5363 pam_tree_clone[i] = new TTree("NeutronD", "PAMELA neutron detector level2 data ");
5364 pam_tree_clone[i]->Branch("NDLevel2", "NDLevel2", GetPointerTo("NDLevel2"));
5365 cout << "NeutronD : branch NDLevel2" << endl;
5366 i++;
5367 };
5368 // Anticounters
5369 if (AC) {
5370 pam_tree_clone[i] = new TTree("Anticounter", "PAMELA anticounter detector level2 data ");
5371 pam_tree_clone[i]->Branch("AcLevel2", "AcLevel2", GetPointerTo("AcLevel2"));
5372 cout << "Anticounter : branch AcLevel2" << endl;
5373 i++;
5374 };
5375 // OrbitalInfo
5376 if (ORB) {
5377 pam_tree_clone[i] = new TTree("OrbitalInfo", "PAMELA orbital info ");
5378 pam_tree_clone[i]->Branch("OrbitalInfo", "OrbitalInfo", GetPointerTo("OrbitalInfo"));
5379 cout << "OrbitalInfo : branch OrbitalInfo" << endl;
5380 if(NUC){
5381 pam_tree_clone[i]->Branch("TrackNuclei","TClonesArray",&orb_nuc_obj);
5382 cout << "OrbitalInfo : branch TrackNuclei" << endl;
5383 }
5384 if(EXT){
5385 pam_tree_clone[i]->Branch("RecoveredTrack","TClonesArray",&orb_ext_obj);
5386 cout << "OrbitalInfo : branch RecoveredTrack" << endl;
5387 if(NUC){
5388 pam_tree_clone[i]->Branch("RecoveredTrackNuclei","TClonesArray",&orb_ext_nuc_obj);
5389 cout << "OrbitalInfo : branch RecoveredTrackNuclei" << endl;
5390 }
5391 }
5392 i++;
5393 };
5394 // GPamela
5395 if (GP) {
5396 pam_tree_clone[i] = new TTree("h20", "GPAMELA info ");
5397 pam_tree_clone[i]->Branch("GPamela", "GPamela", GetPointerTo("GPamela"), 32000, 1);//split
5398 cout << "GPamela : branch GPamela" << endl;
5399 i++;
5400 };
5401
5402
5403 cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" << endl;
5404
5405 }
5406
5407 /**
5408 * Fill tree (created with CreateCloneTrees)
5409 *
5410 */
5411 //void PamLevel2::FillNewPamTree(TTree *T){
5412 void PamLevel2::FillCloneTrees() {
5413
5414 // cout << "PamLevel2::FillCloneTrees()" << irunentry << endl;
5415
5416 for (Int_t i = 0; i < NCLONES; i++) {
5417 if (pam_tree_clone[i])
5418 pam_tree_clone[i]->Fill();
5419 }
5420 if (sel_tree_clone)
5421 sel_tree_clone->Fill();
5422
5423 }
5424
5425 TTree* PamLevel2::GetCloneTree(TString name) {
5426
5427 for (Int_t i = 0; i < NCLONES; i++) {
5428 if (pam_tree_clone[i]) {
5429 TString na = pam_tree_clone[i]->GetName();
5430 if (!name.CompareTo(na))
5431 return pam_tree_clone[i];
5432 };
5433 }
5434 if (run_tree_clone) {
5435 TString na = run_tree_clone->GetName();
5436 if (!name.CompareTo(na))
5437 return run_tree_clone;
5438 }
5439 if (sel_tree_clone) {
5440 TString na = sel_tree_clone->GetName();
5441 if (!name.CompareTo(na))
5442 return sel_tree_clone;
5443 }
5444 if (proc_tree_clone && PROC) {
5445 TString na = proc_tree_clone->GetName();
5446 if (!name.CompareTo(na))
5447 return proc_tree_clone;
5448 }
5449 return NULL;
5450
5451 }
5452 void PamLevel2::WriteCloneTrees() {
5453 cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" << endl;
5454 cout << "Write clones of PAMELA trees " << endl;
5455 cout << run_tree_clone->GetName() << endl;
5456 if (SELLI != 2) {// Emiliano
5457 if (run_tree_clone->GetBranch("DeadLiveTime")->GetEntries() < run_tree->GetEntries())
5458 run_tree_clone->GetBranch("DeadLiveTime")->Fill();
5459 };
5460 run_tree_clone->Write();
5461 if (SELLI != 2) { //Emiliano
5462 cout << sel_tree_clone->GetName() << endl;
5463 sel_tree_clone->Write();
5464 };
5465 for (Int_t i = 0; i < NCLONES; i++) {
5466 if (pam_tree_clone[i]) {
5467 cout << pam_tree_clone[i]->GetName() << endl;
5468 pam_tree_clone[i]->Write(pam_tree_clone[i]->GetName(),TObject::kOverwrite);
5469 };
5470 }
5471
5472 if ( PROC && false ){//EMEMEMEM
5473 proc_tree_clone->Write("ProcessingInfo",TObject::kOverwrite);
5474 }
5475 cout << "+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+" << endl;
5476
5477 }
5478
5479 /**
5480 * Method to get level2-trees entry, the corresponding run entry and (if required) the level0 entry.
5481 */
5482 Int_t PamLevel2::GetEntry(Long64_t iee) {
5483
5484 if (!pam_tree) {
5485 cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- level2 trees not loaded" << endl;
5486 return 0;
5487 }
5488
5489 //
5490 // This is a sort of bug: if you don't have the run tree you don't want to exit here you want to have loaded the event anyway...
5491 //
5492 // if(!run_tree ){
5493 // cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- run tree not loeaded"<<endl;
5494 // return 0;
5495 // }
5496
5497 //-------------------------------
5498 if (!pam_tree->GetEntry(iee)) {
5499 cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- error reading pam tree" << endl;
5500 return 0;
5501 }
5502 //
5503 // ... that's way I put it here. Notice that nothing change in the code (is backward compatible) since in any case you return with 0.
5504 // in theory one would like to return 1 if run is not loaded but now I don't have the will to add that 2 lines of code and it is not
5505 // a problem if you don't check the return code of getentry.
5506 //
5507 if (!RUN || !run_tree) {
5508 if (TRK0 || CAL0 || TOF0 || RUN) { //forse cosi` va bene per tornare 1?
5509 cout << " Int_t PamLevel2::GetEntry(Int_t) -- ERROR -- run tree not loaded" << endl;
5510 return 0;
5511 }
5512 else {
5513 return 1; //cosi` se non c'e` run esce qua...
5514 }
5515 }
5516
5517 //-------------------------------
5518 //
5519 if ( fUpdateRunInfo ) UpdateRunInfo(iee); // Emiliano
5520 if (SELLI == 0 || SELLI == 2) irunentry = iee - runfirstentry;
5521
5522 return 1;
5523
5524 }
5525
5526 TrkLevel0 *PamLevel2::GetTrkLevel0() {
5527 if (!TRK0)
5528 return NULL;
5529 if (!GetYodaEntry()) {
5530 cout << " Int_t PamLevel2::GetTrkLevel0() -- ERROR -- error reading level0 tree" << endl;
5531 return 0;
5532 }
5533 return trk0_obj;
5534 }
5535 ;
5536 CaloLevel0 *PamLevel2::GetCaloLevel0() {
5537 if (!CAL0)
5538 return NULL;
5539 if (!GetYodaEntry()) {
5540 cout << " Int_t PamLevel2::GetCaloLevel0() -- ERROR -- error reading level0 tree" << endl;
5541 return 0;
5542 }
5543 return calo0_obj;
5544 }
5545 ;
5546
5547 /**
5548 * Method to retrieve the level0 tree (YODA tree) that contains the current event.
5549 * Given the run ID (...), if needed it query the DB and load the proper file.
5550 * @return Pointer to the tree
5551 */
5552
5553 TTree* PamLevel2::GetYodaTree() {
5554
5555 // cout << "TTree* PamLevel2::GetYodaTree( )"<<endl;
5556 //===================================
5557 // check if iroot has changed
5558 //===================================
5559 if (irun < 0) {
5560 cout << "PamLevel2::GetYodaTree() -- ERROR -- irun = " << irun << endl;
5561 if (DBG) cout << "In order to use this method you have to first load the RunInfo tree "<<endl;
5562 return NULL;
5563 }
5564 Int_t irootnew = GetRunInfo()->ID_ROOT_L0;
5565 if (DBG){
5566 cout << "iroot "<<iroot<<endl;
5567 cout << "irootnew "<<irootnew<<endl;
5568 }
5569
5570 //===================================
5571 // load the level0 file
5572 // (if not already loaded)
5573 //===================================
5574 if (iroot != irootnew || !l0_tree) {
5575 iroot = irootnew;
5576 //===================================
5577 // open the DB connection
5578 // (if not already opened)
5579 //===================================
5580 if (!dbc || (dbc && !dbc->IsConnected()))
5581 SetDBConnection();
5582 GL_ROOT glroot = GL_ROOT();
5583 if (glroot.Query_GL_ROOT(iroot, dbc)) {
5584 cout << "TTree* PamLevel2::GetYodaTree( ) -- ERROR -- level0 file iroot = " << iroot << " does not exists"
5585 << endl;
5586 return NULL;
5587 };
5588 TString filename = glroot.PATH + glroot.NAME;
5589 if (l0_file) {
5590 l0_file->Close();
5591 l0_file->Delete();
5592 }
5593 cout << "Opening LEVEL0 file: " << filename << endl;
5594 FileStat_t t;
5595 if (gSystem->GetPathInfo(filename.Data(), t)) {
5596 cout << " PamLevel2::GetYodaTree() -- ERROR opening file " << endl;
5597 return NULL;
5598 }
5599 l0_file = new TFile(filename);
5600 if (!l0_file)
5601 return NULL;
5602 l0_tree = (TTree*) l0_file->Get("Physics");
5603 if (!h0_obj)
5604 h0_obj = new EventHeader();
5605 l0_tree->SetBranchAddress("Header", &h0_obj);
5606 yprevshift = 0; // yes, yprevshift is the shift in the level0, prevshift is the shift in the level2
5607 //---------------------------------------------------
5608 // TRACKER:
5609 if (TRK0) {
5610 if (!trk0_obj) {
5611 trk0_obj = new TrkLevel0();
5612 trk0_obj->Set();
5613 };
5614 l0_tree->SetBranchAddress("Tracker", trk0_obj->GetPointerToTrackerEvent());
5615 }
5616 //--------------------------------------------------
5617 // CALORIMETER:
5618 if (CAL0) {
5619 if (!calo0_obj) {
5620 calo0_obj = new CaloLevel0();
5621 calo0_obj->Set();
5622 };
5623 l0_tree->SetBranchAddress("Calorimeter", calo0_obj->GetPointerToCalorimeterEvent());
5624 }
5625 //---------------------------------------------------
5626 // TOF:
5627 if (TOF0) {
5628 cout << "PamLevel2::GetYodaTree() --- level0 TOF not implemented " << endl;
5629 }
5630
5631 dbc->Close(); // EMILIANO, do not leave open connections, open only when needed
5632 delete dbc;
5633 dbc=0;
5634
5635 };
5636
5637 if (TRK0) {
5638 if(!dbc || (dbc && !dbc->IsConnected()))SetDBConnection();
5639 TrkParams::SetCalib(run_obj, dbc);
5640 TrkParams::LoadCalib();
5641 if (!TrkParams::CalibIsLoaded()) {
5642 cout << " TTree* PamLevel2::GetYodaTree( ) -- WARNING -- Calibration not loaded" << endl;
5643 };
5644 if(dbc){
5645 dbc->Close(); // EMILIANO, do not leave open connections, open only when needed
5646 delete dbc;
5647 dbc=0;
5648 };
5649 }
5650 return l0_tree;
5651 }
5652
5653 /**
5654 * Method to retrieve the level0 tree (YODA tree) that contains the current event.
5655 */
5656 Int_t PamLevel2::GetYodaEntry() {
5657
5658 Long64_t iev = this->GetReadEntry();
5659
5660 // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
5661 // if it is a full file (not preselected)
5662 // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
5663 // if (SELLI == 0 || SELLI == 2 || !hasL0EE) {
5664
5665 if (!GetYodaTree()){
5666 printf(" PamLevel2::GetYodaEntry() : ERROR no level0 file loaded!\n");
5667 return 0;
5668 }
5669
5670 if (irunentry < 0) {
5671 if (DBG) cout << "Int_t PamLevel2::GetYodaEntry() -- ATTENZIONE -- irunentry negativo?!?! "<<(Int_t)irunentry<<endl;
5672 // cout << "Int_t PamLevel2::GetYodaEntry() -- ATTENZIONE -- irunentry negativo?!?! "<<(Int_t)irunentry<<endl; // TOGLITOGLI
5673 irunentry = 0LL;
5674 }
5675 // ---------------------------------
5676 // if file is NOT a preselected file
5677 // ---------------------------------
5678 Long64_t quellagiusta = irunentry + (Long64_t)(run_obj->EV_FROM); // prevshift already included in irunentry
5679
5680 if (DBG){
5681 cout << " irun "<< irun << " irunentry "<< irunentry<<" run_obj->EV_FROM "<<run_obj->EV_FROM <<" quella giusta "<<quellagiusta << endl;
5682 cout << " iroot "<<iroot<<" run_obj->ID_ROOT_L0 "<<run_obj->ID_ROOT_L0<<endl;
5683 cout << " time "<< abstime << endl;
5684 }
5685
5686 ULong64_t obt = 0;
5687 ULong64_t pktn = 0;
5688 if (GetOrbitalInfo()) {
5689 obt = GetOrbitalInfo()->OBT;
5690 pktn = GetOrbitalInfo()->pkt_num;
5691 }
5692
5693 if (!GetOrbitalInfo() && !ISGP) {
5694 cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- missing OrbitalInfo " << endl;
5695 return 0;
5696 }
5697 if (obt == 0 && pktn == 0 && !ISGP) {
5698 cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- level2 event corrupted ?? " << endl;
5699 return 0;
5700 }
5701
5702 // ---------------------------------------------------------------------
5703 // ATTENTION!!!
5704 // If data are not pre-processed with cleaner, the level0 tree may contain
5705 // spurious nested physics packets.
5706 // The GL_RUN variables EV_FROM, EV_TO and NEVENTS counts also these entries
5707 // while level2 tree DOES NOT!!
5708 // This means that "quellagiusta" in these cases is not correct.
5709 // In order to retrieve the correct level0 event, I implemented a check
5710 // of the OBT and pkt-number. In case of mismatch, the level0 entry number
5711 // is shift forward until when the packets match.
5712 // ---------------------------------------------------------------------
5713 Long64_t shift = 0LL;
5714 Long64_t answer = quellagiusta + shift + yprevshift;
5715 Int_t readl0 = 0;
5716 readl0 = l0_tree->GetEntry(answer); // prevshift already included in irunentry
5717
5718 if (DBG){
5719 printf(" siamo qui shift %lld yprevshift %lld answer %lld \n",shift,yprevshift,answer);
5720 }
5721
5722
5723 if (ISGP) {
5724 obt = (UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()); //BARBATRUCCO
5725 pktn = (UInt_t)(GetEventHeader()->GetPscuHeader()->GetCounter()); //BARBATRUCCO
5726 }
5727
5728 while ( (obt != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()) || pktn != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetCounter())) && (quellagiusta + (Long64_t) shift) < GetYodaTree()->GetEntries() && shift < maxshift ){
5729 if ( isSync && shift == 0LL ){
5730 printf(" PamLevel2::GetYodaEntry() ERROR! sync file but the level0 entry not found in place!!! \n");
5731 cout << " OK this is a bug, write to Emiliano, Emiliano.Mocchiutti@ts.infn.it " << endl;
5732 cout << "\nFor bug reporting instructions, please see for example:\n";
5733 cout << " <http://www.ts.infn.it/~mocchiut/bugs/bugs.html>.\n";
5734 cout << " " << endl;
5735 }
5736 if (shift > 0) {
5737 if (DBG) cout << " PKTNUM L2 --- " << pktn << " --- L0 --- " << GetEventHeader()->GetPscuHeader()->GetCounter() << endl;
5738 if (DBG)
5739 cout << " RUN: ID " << GetRunInfo()->ID << " ID_ROOT_L0 " << run_obj->ID_ROOT_L0 << " ID_RUN_FRAG "
5740 << GetRunInfo()->ID_RUN_FRAG << " EV_FROM " << GetRunInfo()->EV_FROM << endl;
5741 if (DBG)
5742 cout << " L2 <--> L0 mismatch ( irun " << irun << " irunentry " << irunentry << " shift " << shift
5743 << " prevshift " << prevshift << " )" << endl;
5744 }
5745 answer = quellagiusta + shift+ yprevshift;
5746 readl0 = l0_tree->GetEntry(answer);
5747 // printf(" inside while shift %lld yprevshift %lld answer %lld \n",shift,yprevshift,answer);//TOGLITOGLI
5748
5749 if (!GetEventHeader()) {
5750 cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- missing EventHeader " << endl;
5751 return 0;
5752 }
5753
5754 //
5755 if (yprevshift != 0 && (quellagiusta + (Long64_t) shift) == GetYodaTree()->GetEntries()) {
5756 if (DBG) printf(" reset inside while shift %lld yprevshift %lld answer %lld \n",shift,yprevshift,answer);
5757 // printf(" reset inside while shift %lld yprevshift %lld answer %lld \n",shift,yprevshift,answer);//TOGLITOGLI
5758 yprevshift = 0LL;
5759 shift = -1LL;
5760 };
5761
5762 shift++;
5763 }
5764
5765
5766 if ( obt != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()) || pktn != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetCounter()) ){
5767 if ( isSync ){
5768 printf(" PamLevel2::GetYodaEntry() ERROR! sync file but the level0 entry not found AT ALL!!! \n");
5769 cout << " OK this is a bug, write to Emiliano, Emiliano.Mocchiutti@ts.infn.it " << endl;
5770 cout << "\nFor bug reporting instructions, please see for example:\n";
5771 cout << " <http://www.ts.infn.it/~mocchiut/bugs/bugs.html>.\n";
5772 cout << " " << endl;
5773 }
5774 cout << "Int_t PamLevel2::GetYodaEntry() -- WARNING -- " << endl;
5775 cout << " Big trouble here, no such event in Level0 data! (NB maxshift set to " << maxshift << " )" << endl;
5776 cout << " Nested and/or DarthVader skipped packets in fragmented run? checking and trying to fix " <<endl;
5777 // query the DB for runs containing the event, loop over LEVEL0 files which could contain the level0 event and try to find it
5778 // ma nel mezzo del cammin per ogni run che pesco devo vedere la posizione relativa di iev rispetto al runheader nel livello2 per andare a cercare nel posto giusto
5779 // connect to db
5780 if (!dbc || (dbc && !dbc->IsConnected())) SetDBConnection(); //Emiliano
5781 //
5782 if (GetOrbitalInfo()){
5783 abstime = GetOrbitalInfo()->absTime;
5784 } else {
5785 printf(" PamLevel2::GetYodaEntry() ERROR! no OrbitalInfo, cannot get the absolute time for event \n");
5786 return 0;
5787 }
5788 // query DB looking for runs containing the processed event
5789 TSQLResult *pResult;
5790 TSQLRow *Row = NULL;
5791 TString myquery = Form("select ID,NEVENTS from GL_RUN where RUNHEADER_TIME<=%lld and RUNTRAILER_TIME>=%lld;",abstime,abstime);
5792 if ( DBG ) printf(" query is %s \n",myquery.Data());
5793 // printf(" query is %s \n",myquery.Data());// TOGLITOGLI
5794 pResult = dbc->Query(myquery.Data());
5795 if (!pResult->GetRowCount()){
5796 printf(" PamLevel2::GetYodaEntry() ERROR! event is not included in any run!!! \n");
5797 cout << " OK this is a bug, write to Emiliano, Emiliano.Mocchiutti@ts.infn.it " << endl;
5798 cout << "\nFor bug reporting instructions, please see for example:\n";
5799 cout << " <http://www.ts.infn.it/~mocchiut/bugs/bugs.html>.\n";
5800 cout << " " << endl;
5801 return 0;
5802 }
5803 if ( pResult->GetRowCount() == 1 ){
5804 if (DBG) printf(" PamLevel2::GetYodaEntry() - WARNING - YodaEntry not found but only one run containing the event, it should not happen \n");
5805 // printf(" PamLevel2::GetYodaEntry() - WARNING - YodaEntry not found but only one run containing the event, it should not happen \n");//TOGLITOGLI
5806 }
5807 for( Int_t ru=0; ru < pResult->GetRowCount(); ru++){ // loop over runs containing the event
5808 if (Row) delete Row;
5809 Row = pResult->Next();
5810 if( Row == NULL ) break;
5811 UInt_t idrun = (UInt_t)atoll(Row->GetField(0));
5812 UInt_t nev = (UInt_t)atoll(Row->GetField(1));
5813 if (DBG) printf(" inside loop over runs: ru %i idrun %i nev %i \n",ru,idrun,nev);
5814 // printf(" inside loop over runs: ru %i idrun %i nev %i \n",ru,idrun,nev);//TOGLITOGLI
5815
5816 // now look for this run in the level2 file, it must be present! code is taken from updateruninfo of course
5817 Bool_t rfound = false;
5818 totrunentry = 0LL;
5819 runfirstentry = 0LL;
5820 for (Int_t r=0; r< run_tree->GetEntries();r++){
5821 run_tree->GetEntry(r);//update runinfo
5822 if ( r > 0 ){
5823 totrunentrymin = totrunentrymax+1;
5824 } else {
5825 totrunentrymin = 0LL;
5826 }
5827 totrunentry += GetRunInfo()->NEVENTS;
5828 totrunentrymax = totrunentry - 1 - prevshift; // prevshift is needed to handle nested+DV_skipped packets
5829 irun = r;
5830
5831 if (idrun == GetRunInfo()->ID){
5832 if ( totrunentrymin > iev ){ // there is a shift (nested+DV_skipped packets)
5833 if (DBG) printf("PamLevel2::GetYodaEntry - unconsistent iev - nevents, probable DBL0-L2 async\n");
5834 if (DBG) printf("PamLevel2::GetYodaEntry - totrunentrymin %lld iev %lld prevshift %lld totrunentrymax %lld \n",totrunentrymin,iev,prevshift,totrunentrymax);
5835 // printf("PamLevel2::GetYodaEntry - unconsistent iev - nevents, probable DBL0-L2 async\n");
5836 // printf("PamLevel2::GetYodaEntry - totrunentrymin %lld iev %lld prevshift %lld totrunentrymax %lld \n",totrunentrymin,iev,prevshift,totrunentrymax);//TOGLITOGLI
5837 prevshift += (totrunentrymin-iev); // add the new shift to total shift
5838 totrunentrymin -= (totrunentrymin-iev); // shift run position min
5839 totrunentrymax -= (totrunentrymin-iev); // shift run position max
5840 if (DBG) printf("PamLevel2::GetYodaEntry - totrunentrymin %lld iev %lld prevshift %lld totrunentrymax %lld \n",totrunentrymin,iev,prevshift,totrunentrymax);
5841 // printf("PamLevel2::GetYodaEntry - totrunentrymin %lld iev %lld prevshift %lld totrunentrymax %lld \n",totrunentrymin,iev,prevshift,totrunentrymax);//TOGLITOGLI
5842 }
5843 runfirstentry = totrunentrymin; // first entry of the run in the level2
5844
5845
5846 //
5847 if (gltsync)
5848 delete gltsync; // Emiliano
5849 if (!dbc || (dbc && !dbc->IsConnected()))
5850 SetDBConnection(); //Emiliano
5851 gltsync = new GL_TIMESYNC(GetRunInfo()->ID_ROOT_L0, "ID", dbc, false); // Emiliano
5852 if (dbc){
5853 dbc->Close(); // Emiliano
5854 delete dbc;
5855 dbc=0;
5856 }
5857 if (gltsync->DBobt(GetRunInfo()->RUNHEADER_OBT) > gltsync->DBobt(GetRunInfo()->RUNTRAILER_OBT)) { // Emiliano
5858 cout << "Bool_t PamLevel2::UpdateRunInfo(Long64_t iev) -- WARNING -- irun " << irun
5859 << " has RUNHEADER_OBT>=RUNTRAILER_OBT " << endl;
5860 cout
5861 << " (NB!! in this case some events could be assigned to a wrong run)"
5862 << endl;
5863 }
5864 //
5865 if (DBG) printf(" found \n");
5866 // printf(" found \n");//TOGLITOGLI
5867 rfound = true;
5868 //
5869 break;
5870 }
5871 } // loop over run
5872 if ( !rfound ){
5873 printf(" PamLevel2::GetYodaEntry() ERROR! run is not present in the level2 file!!! \n");
5874 cout << " OK this is a bug, write to Emiliano, Emiliano.Mocchiutti@ts.infn.it " << endl;
5875 cout << "\nFor bug reporting instructions, please see for example:\n";
5876 cout << " <http://www.ts.infn.it/~mocchiut/bugs/bugs.html>.\n";
5877 cout << " " << endl;
5878 return 0;
5879 }
5880
5881 // here we got the first run and we can check if it contains the level0 event
5882 if (!GetYodaTree()){
5883 printf(" PamLevel2::GetYodaEntry() : ERROR no level0 file loaded!\n");
5884 return 0;
5885 }
5886
5887 // get the current run entry
5888 irunentry = iev - runfirstentry;
5889 if (irunentry < 0) {
5890 if (DBG) cout << "Int_t PamLevel2::GetYodaEntry() -- ATTENZIONE -- irunentry negativo?!?! "<<(Int_t)irunentry<<endl;
5891 // cout << "Int_t PamLevel2::GetYodaEntry() -- ATTENZIONE -- irunentry negativo?!?! "<<(Int_t)irunentry<<endl; // TOGLITOGLI
5892 irunentry = 0LL;
5893 }
5894 // ---------------------------------
5895 // if file is NOT a preselected file
5896 // ---------------------------------
5897 quellagiusta = irunentry + (Long64_t)(run_obj->EV_FROM); // prevshift already included in irunentry
5898
5899 if (DBG){
5900 cout << " irun "<< irun << " irunentry "<< irunentry<<" run_obj->EV_FROM "<<run_obj->EV_FROM <<" quella giusta "<<quellagiusta << endl;
5901 cout << " iroot "<<iroot<<" run_obj->ID_ROOT_L0 "<<run_obj->ID_ROOT_L0<<endl;
5902 cout << " time "<< abstime << endl;
5903 }
5904 // cout << " irun "<< irun << " irunentry "<< irunentry<<" run_obj->EV_FROM "<<run_obj->EV_FROM <<" quella giusta "<<quellagiusta << endl;
5905 // cout << " iroot "<<iroot<<" run_obj->ID_ROOT_L0 "<<run_obj->ID_ROOT_L0<<endl;
5906 // cout << " time "<< abstime << endl; // TOGLITOGLI
5907
5908 shift = 0;
5909 answer = quellagiusta + shift + yprevshift; // prevshift already included in irunentry
5910 readl0 = l0_tree->GetEntry(answer); // prevshift already included in irunentry
5911
5912 if (DBG){
5913 printf(" siamo qua shift %lld yprevshift %lld answer %lld \n",shift,yprevshift,answer);
5914 }
5915 // printf(" siamo qua shift %lld yprevshift %lld answer %lld \n",shift,yprevshift,answer);//TOGLITOGLI
5916
5917 while ( (obt != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()) || pktn != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetCounter())) && (quellagiusta + (Long64_t) shift) < GetYodaTree()->GetEntries() && shift < maxshift ){
5918 if (shift > 0) {
5919 if (DBG) cout << " PKTNUM L2 --- " << pktn << " --- L0 --- " << GetEventHeader()->GetPscuHeader()->GetCounter() << endl;
5920 if (DBG)
5921 cout << " RUN: ID " << GetRunInfo()->ID << " ID_ROOT_L0 " << run_obj->ID_ROOT_L0 << " ID_RUN_FRAG "
5922 << GetRunInfo()->ID_RUN_FRAG << " EV_FROM " << GetRunInfo()->EV_FROM << endl;
5923 if (DBG)
5924 cout << " L2 <--> L0 mismatch ( irun " << irun << " irunentry " << irunentry << " shift " << shift
5925 << " prevshift " << prevshift << " )" << endl;
5926 }
5927 answer = quellagiusta + shift+ yprevshift;
5928 readl0 = l0_tree->GetEntry(answer);
5929 // printf(" inside inside while shift %lld yprevshift %lld answer %lld \n",shift,yprevshift,answer);//TOGLITOGLI
5930
5931 if (!GetEventHeader()) {
5932 cout << "Int_t PamLevel2::GetYodaEntry() -- ERROR -- missing EventHeader " << endl;
5933 return 0;
5934 }
5935 //
5936 if (yprevshift != 0 && (quellagiusta + (Long64_t) shift) == GetYodaTree()->GetEntries()) {
5937 if (DBG) printf(" reset inside while shift %lld yprevshift %lld answer %lld \n",shift,yprevshift,answer);
5938 // printf(" reset inside while shift %lld yprevshift %lld answer %lld \n",shift,yprevshift,answer);//TOGLITOGLI
5939 yprevshift = 0;
5940 shift = -1;
5941 };
5942
5943 shift++;
5944 }
5945
5946 if ( obt != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()) || pktn != (UInt_t)(GetEventHeader()->GetPscuHeader()->GetCounter()) ){
5947 //still not the good run... continue with the nex one!
5948 printf("still not the good run... continue with the nex one!\n");
5949 } else {
5950 if (DBG) cout << "LA ENTRY GIUSTA E`: "<<quellagiusta<<" (spero...)"<<endl;
5951 // cout << "LA ENTRY GIUSTA E`: "<<answer<<" (spero...)"<<endl;//TOGLITOGLI
5952 if (shift > 1) yprevshift = (shift - 1);
5953 if (Row) delete Row;
5954 delete pResult;
5955 if (dbc){
5956 dbc->Close(); // Emiliano
5957 delete dbc;
5958 dbc=0;
5959 }
5960 il0entry = answer;
5961 return readl0;
5962 }
5963 // perhaps it is all
5964 }// loop over runs containing the event
5965 if (Row) delete Row;
5966 delete pResult;
5967 if (dbc){
5968 dbc->Close(); // Emiliano
5969 delete dbc;
5970 dbc=0;
5971 }
5972 // arriving here it means no run found, cannot be! error!
5973 printf(" PamLevel2::GetYodaEntry() ERROR! run is not present in the level0 files!!! \n");
5974 cout << " OK this is a bug, write to Emiliano, Emiliano.Mocchiutti@ts.infn.it " << endl;
5975 cout << "\nFor bug reporting instructions, please see for example:\n";
5976 cout << " <http://www.ts.infn.it/~mocchiut/bugs/bugs.html>.\n";
5977 cout << " " << endl;
5978 return 0;
5979 } else {
5980 if (DBG) cout << "=> LA ENTRY GIUSTA E`: "<<answer<<" (spero...)"<<endl;
5981 // cout << "=> LA ENTRY GIUSTA E`: "<<answer<<" (spero...)"<<endl;
5982 // printf("obt %lld (UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()) %i pktn %lld (UInt_t)(GetEventHeader()->GetPscuHeader()->GetCounter() %i \n",obt,(UInt_t)(GetEventHeader()->GetPscuHeader()->GetOrbitalTime()), pktn, (UInt_t)(GetEventHeader()->GetPscuHeader()->GetCounter()) );
5983 if (shift > 1) yprevshift = (shift - 1);
5984 il0entry = answer;
5985 return readl0;
5986 }
5987
5988 /* } // if selli 0 || 2
5989 if ( SELLI == 1 && hasL0EE ){
5990 sel_tree->GetEntry(iev);
5991 Long64_t answer = il0entry;
5992 Int_t readl0 = 0;
5993 readl0 = l0_tree->GetEntry(answer);
5994 return readl0;
5995 }*/
5996 printf(" PamLevel2::GetYodaEntry() ERROR! \n");
5997 cout << " Entry not found! OK this is a bug, write to Emiliano, Emiliano.Mocchiutti@ts.infn.it " << endl;
5998 cout << "\nFor bug reporting instructions, please see for example:\n";
5999 cout << " <http://www.ts.infn.it/~mocchiut/bugs/bugs.html>.\n";
6000 cout << " " << endl;
6001 return 0;
6002 }
6003
6004 /**
6005 * \Brief Set DB connection
6006 */
6007 Bool_t PamLevel2::SetDBConnection() {
6008
6009 // cout << "PamLevel2::SetDBConnection()" << endl;
6010 if (DBG) {
6011 cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << endl;
6012 cout << "Connecting to DB" << endl;
6013 cout << "HOST " << host << endl;
6014 cout << "USER " << user << endl;
6015 cout << "PSW " << psw << endl;
6016 }
6017 Bool_t notconn = true;
6018 Int_t trials = 0;
6019 while ( notconn && trials < 10 ){
6020 // gSystem->Sleep(500);
6021 dbc = TSQLServer::Connect(host.Data(), user.Data(), psw.Data());
6022 //dbc->Connect(host.Data(), user.Data(), psw.Data());
6023 if ( dbc ) notconn = false;
6024 if (DBG) printf("<%i> test connection...\n ",trials);
6025 if (!dbc){
6026 if (DBG) printf(" :( failed, no pointer \n");
6027 notconn = true;
6028 // return false;
6029 };
6030 if (dbc && !dbc->IsConnected()){
6031 if (DBG) printf(" :( failed, no connection \n");
6032 notconn = true;
6033 // return false;
6034 };
6035 trials++;
6036 };
6037 if ( notconn ) return false;
6038 //
6039 if (DBG) printf("=connected!\n");
6040 stringstream myquery; // EMILIANO
6041 myquery.str(""); // EMILIANO
6042 myquery << "SET time_zone='+0:00'"; // EMILIANO
6043 dbc->Query(myquery.str().c_str()); // EMILIANO
6044 if ( DBG ) printf("~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n");
6045 return true;
6046
6047 }
6048
6049 /**
6050 * \Brief Add a friend to the pamela chain.
6051 * @param cname name of the chain to be added
6052 */
6053
6054 TChain* PamLevel2::AddFriend(TString cname) {
6055
6056 if (!GetPamTree()) {
6057 cout << " TChain* PamLevel2::AddFriend(TString cname) --- a pamela tree must be created first" << endl;
6058 return NULL;
6059 }
6060
6061 TChain *c = new TChain(cname.Data());
6062
6063 TIter next(GetPamTree()->GetListOfFiles());
6064 Int_t nf = 0;
6065 TChainElement* element = 0;
6066 while ((element = (TChainElement*) next())) {
6067 c->Add(element->GetTitle());
6068 nf++;
6069 }
6070
6071 GetPamTree()->AddFriend(cname.Data());
6072
6073 cout << "external chain created and added to pamela friends :" << cname << endl;
6074 cout << "n.files " << nf << endl;
6075
6076 return c;
6077
6078 }
6079
6080 /**
6081 * Returns the current read entry. This method simply returns the result of the call to
6082 * pam_tree->GetReadEntry(), so it is entirely handled by ROOT.
6083 */
6084 Long64_t PamLevel2::GetReadEntry() {
6085 return pam_tree->GetReadEntry();
6086 }
6087
6088 /**
6089 * Sets the sorting method. If the new method is different from the previous, the issorted
6090 * flag is set to false, forcing a new sort the next time GetTrack is called.
6091 * @see GetTrack
6092 */
6093 void PamLevel2::SetSortingMethod(TString how) {
6094 if (howtosort != how) {
6095 issorted = false;
6096 issorted_new = false;
6097 }
6098 howtosort = how;
6099 }
  ViewVC Help
Powered by ViewVC 1.1.23