| 1 |
|
| 2 |
#if !defined(__CINT__) || defined(__MAKECINT__) |
| 3 |
|
| 4 |
|
| 5 |
// ROOT includes // |
| 6 |
#include <TString.h> |
| 7 |
#include <TH1F.h> |
| 8 |
#include <TH2F.h> |
| 9 |
|
| 10 |
// standard libraries includes // |
| 11 |
#include <stdlib.h> |
| 12 |
#include <iostream> |
| 13 |
//#include <fstream> |
| 14 |
#include <iomanip> |
| 15 |
using namespace std; |
| 16 |
|
| 17 |
// PAMELA includes // |
| 18 |
#include <PamLevel2.h> |
| 19 |
|
| 20 |
#endif |
| 21 |
|
| 22 |
//====================== |
| 23 |
// global variables |
| 24 |
//====================== |
| 25 |
Int_t tot=0; |
| 26 |
Int_t sel[]={0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}; |
| 27 |
|
| 28 |
TrkLevel2 *trkl2 = 0; |
| 29 |
CaloLevel2 *calol2 = 0; |
| 30 |
ToFLevel2 *tofl2 = 0; |
| 31 |
TrigLevel2 *trigl2 = 0; |
| 32 |
AcLevel2 *acl2 = 0; |
| 33 |
|
| 34 |
//====================== |
| 35 |
// HISTOGRAMS declaration |
| 36 |
//====================== |
| 37 |
|
| 38 |
TH2F* hxytof[6]; |
| 39 |
|
| 40 |
void Dump(ToFTrkVar* tof){ |
| 41 |
|
| 42 |
cout << endl<< "n. tdc hits :"<<tof->npmttdc; |
| 43 |
cout << endl<< "PMTs :"; |
| 44 |
for(Int_t ip=0; ip< tof->npmttdc; ip++)cout << " "<< tof->pmttdc[ip]; |
| 45 |
cout << endl<< "n. adc hits :"<<tof->npmtadc; |
| 46 |
cout << endl<< "PMTs :"; |
| 47 |
for(Int_t ip=0; ip< tof->npmtadc; ip++)cout << " "<< tof->pmtadc[ip]; |
| 48 |
cout << endl<< "dedx :"; |
| 49 |
for(Int_t ib=0; ib< tof->npmtadc; ib++)cout << " "<< tof->dedx[ib]; |
| 50 |
cout << endl<< "beta :"; |
| 51 |
for(Int_t ib=0; ib<13 ; ib++)cout << " "<< tof->beta[ib]; |
| 52 |
cout << endl<< "xtofpos :"; |
| 53 |
for(Int_t ib=0; ib<3 ; ib++)cout << " "<< tof->xtofpos[ib]; |
| 54 |
cout << endl<< "ytofpos :"; |
| 55 |
for(Int_t ib=0; ib<3 ; ib++)cout << " "<<tof->ytofpos[ib]; |
| 56 |
cout <<endl<<"============================================="<<endl; |
| 57 |
|
| 58 |
} |
| 59 |
|
| 60 |
|
| 61 |
void Dump(ToFLevel2* l2){ |
| 62 |
cout<<"ToF event: "<<l2->npmt()<<" hit PMTs"<<endl; |
| 63 |
cout<<" id name plane paddle pmt adc tdc"<<endl; |
| 64 |
for(int ip=0; ip<l2->npmt(); ip++){ |
| 65 |
Int_t iplane = -1; |
| 66 |
Int_t ipaddle = -1; |
| 67 |
Int_t ipmt = -1; |
| 68 |
// cout << l2->GetToFPMT(ip) << endl; |
| 69 |
Int_t id = l2->GetToFPMT(ip)->pmt_id; |
| 70 |
float adc = l2->GetToFPMT(ip)->adc; |
| 71 |
float tdc = l2->GetToFPMT(ip)->tdc_tw; // temporaneo |
| 72 |
cout <<setw(3)<<id<< setw(10) << l2->GetPMTName(id,iplane,ipaddle,ipmt) << setw(7) << iplane << setw(7)<< ipaddle << setw(7) << ipmt << setw(7) << adc << setw(7) << tdc <<endl; |
| 73 |
} |
| 74 |
cout <<"============================================="<<endl; |
| 75 |
|
| 76 |
} |
| 77 |
|
| 78 |
void Dump(TrigLevel2* l2){ |
| 79 |
cout << "Trigger Configuration (hex): " << hex; |
| 80 |
cout << l2->trigconf << endl; |
| 81 |
cout <<"============================================="<<endl; |
| 82 |
cout << "Trigger Pattern (hex): " << hex; |
| 83 |
for(int i=0; i<6; i++)cout << setw(5) << l2->patterntrig[i]; |
| 84 |
cout << dec << endl; |
| 85 |
cout <<"============================================="<<endl; |
| 86 |
} |
| 87 |
|
| 88 |
void EvaluatePatternTrig(UInt_t *patterntrig, ToFLevel2* l2){ |
| 89 |
|
| 90 |
// TString photoS[48] = { |
| 91 |
// "S11_1A", "S11_1B", "S11_2A", "S11_2B", "S11_3A", "S11_3B", "S11_4A", |
| 92 |
// "S11_4B", |
| 93 |
// "S11_5A", "S11_5B", "S11_6A", "S11_6B", "S11_7A", "S11_7B", "S11_8A", |
| 94 |
// "S11_8B", |
| 95 |
// "S12_1A", "S12_1B", "S12_2A", "S12_2B", "S12_3A", "S12_3B", "S12_4A", |
| 96 |
// "S12_4B", "S12_5A", "S12_5B", "S12_6A", "S12_6B", |
| 97 |
// "S21_1A", "S21_1B", "S21_2A", "S21_2B", |
| 98 |
// "S22_1A", "S22_1B", "S22_2A", "S22_2B", |
| 99 |
// "S31_1A", "S31_1B", "S31_2A", "S31_2B", "S31_3A", "S31_3B", |
| 100 |
// "S32_1A", "S32_1B", "S32_2A", "S32_2B", "S32_3A", "S32_3B" |
| 101 |
// }; |
| 102 |
// patterntrig(3) --> S3 [49 to 38] 12 bits: (S32b,S32a,S31b,S31a) |
| 103 |
// patterntrig(4) --> S2 [37 to 30] 8 bits: (S22b,S22a,S21b,S21a) |
| 104 |
// patterntrig(5) --> S12 [29 to 18] 12 bits: (S12b,S12a) |
| 105 |
// patterntrig(6) --> S11 [17 to 2] 16 bits: (S11b,S11a) |
| 106 |
|
| 107 |
// UInt_t patterntrig[6]; |
| 108 |
for(int i=0; i<6; i++)*(patterntrig+i)=0; |
| 109 |
|
| 110 |
UInt_t word[]={ |
| 111 |
5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5, //S11 |
| 112 |
4,4,4,4,4,4,4,4,4,4,4,4, //S12 |
| 113 |
3,3,3,3,3,3,3,3, //S2 |
| 114 |
2,2,2,2,2,2,2,2,2,2,2,2 //S3 |
| 115 |
}; |
| 116 |
UInt_t bit[] ={ |
| 117 |
0,8,1,9,2,10,3,11,4,12,5,13,6,14,7,15, //S11 |
| 118 |
0,6,1,7,2,8,3,9,4,10,5,11, //S12 |
| 119 |
0,2,1,3,4,6,5,7, //S2 |
| 120 |
0,3,1,4,2,5,6,9,7,10,8,11 //S3 |
| 121 |
}; |
| 122 |
|
| 123 |
for(int ip=0; ip<l2->npmt(); ip++){ |
| 124 |
Int_t iplane = -1; |
| 125 |
Int_t ipaddle = -1; |
| 126 |
Int_t ipmt = -1; |
| 127 |
// cout << l2->GetToFPMT(ip) << endl; |
| 128 |
Int_t id = l2->GetToFPMT(ip)->pmt_id; |
| 129 |
float tdc = l2->GetToFPMT(ip)->tdc_tw; // temporaneo |
| 130 |
|
| 131 |
if(tdc<4095) *(patterntrig+word[id]) = *(patterntrig+word[id]) | (UInt_t)(pow(2.,(float)bit[id])); |
| 132 |
|
| 133 |
} |
| 134 |
|
| 135 |
|
| 136 |
|
| 137 |
|
| 138 |
} |
| 139 |
|
| 140 |
//=============================================================================== |
| 141 |
// |
| 142 |
// |
| 143 |
// |
| 144 |
// |
| 145 |
//=============================================================================== |
| 146 |
bool Select( PamLevel2* event ){ |
| 147 |
// return false; |
| 148 |
TrkParams::SetVerboseMode(); |
| 149 |
|
| 150 |
//------------------------------------------- |
| 151 |
// check if the event contains selection info |
| 152 |
//------------------------------------------- |
| 153 |
if( event->GetTrkLevel2()==0 ){ |
| 154 |
cout << "Missing TrkLevel2 object "<<endl; |
| 155 |
return false; |
| 156 |
} |
| 157 |
if( event->GetToFLevel2()==0 ){ |
| 158 |
cout << "Missing ToFLevel2 object "<<endl; |
| 159 |
return false; |
| 160 |
} |
| 161 |
if( event->GetCaloLevel2()==0 ){ |
| 162 |
cout << "Missing CaloLevel2 object "<<endl; |
| 163 |
return false; |
| 164 |
} |
| 165 |
if( event->GetTrigLevel2()==0 ){ |
| 166 |
cout << "Missing TrigLevel2 object "<<endl; |
| 167 |
return false; |
| 168 |
} |
| 169 |
if( event->GetAcLevel2()==0 ){ |
| 170 |
cout << "Missing AcLevel2 object "<<endl; |
| 171 |
return false; |
| 172 |
} |
| 173 |
if( event->GetRunInfo()==0 ){ |
| 174 |
cout << "Missing RunInfo object "<<endl; |
| 175 |
return false; |
| 176 |
} |
| 177 |
|
| 178 |
trkl2 = event->GetTrkLevel2(); |
| 179 |
calol2 = event->GetCaloLevel2(); |
| 180 |
tofl2 = event->GetToFLevel2(); |
| 181 |
trigl2 = event->GetTrigLevel2(); |
| 182 |
acl2 = event->GetAcLevel2(); |
| 183 |
|
| 184 |
tot++; |
| 185 |
|
| 186 |
|
| 187 |
//------------------------------------------------------------------ |
| 188 |
// Orbital selection |
| 189 |
//------------------------------------------------------------------ |
| 190 |
bool ORB__OK = false; |
| 191 |
if( |
| 192 |
event->GetOrbitalInfo()->Babs>0.22 && //outside SSA |
| 193 |
true ) ORB__OK = true; |
| 194 |
if( !ORB__OK )return false; |
| 195 |
sel[0]++; |
| 196 |
|
| 197 |
//--------------------------------------------------------- |
| 198 |
// trigger configuration: TOF 1 (standard) |
| 199 |
//--------------------------------------------------------- |
| 200 |
bool trigg_patt_ok = false; |
| 201 |
if ( trigl2->trigconf & (1<<0) ) trigg_patt_ok = true; |
| 202 |
if(!trigg_patt_ok)return false; |
| 203 |
sel[1]++; |
| 204 |
|
| 205 |
//--------------------------------------------------------- |
| 206 |
// single track |
| 207 |
//--------------------------------------------------------- |
| 208 |
if( trkl2->GetNTracks()!=1 ) return false; // singola traccia (fisica) |
| 209 |
event->SetSortingMethod("+CAL"); // usa solo l'info del calorimetro per scegliere l'immagine |
| 210 |
PamTrack *track = event->GetTrack(0);// qui viene scelta l'"immagine" fisica |
| 211 |
// if( track->GetScore() <=0 )cout << " >>> Strange track <<< P-score "<< track->GetPScore() << " I-score "<<track->GetIScore()<<endl; |
| 212 |
|
| 213 |
if(!track)return false; |
| 214 |
|
| 215 |
//---------------------------------------------------------- |
| 216 |
// (CHI**2 standard) |
| 217 |
//---------------------------------------------------------- |
| 218 |
float p0 = 1.111588e+00; |
| 219 |
float p1 = 1.707656e+00; |
| 220 |
float p2 = 1.489693e-01; |
| 221 |
float def_0 = 0.07; |
| 222 |
float chi2m025_0 = p0 + fabs(def_0)*p1 + def_0*def_0*p2; |
| 223 |
|
| 224 |
//------------------------------------------------------------------ |
| 225 |
// tracker pre-selection |
| 226 |
//------------------------------------------------------------------ |
| 227 |
TrkTrack *trk = track->GetTrkTrack(); |
| 228 |
|
| 229 |
float chi2m025 = p0 + fabs(trk->GetDeflection())*p1 + trk->GetDeflection()*trk->GetDeflection()*p2; |
| 230 |
float chi2m = pow( chi2m025-chi2m025_0+pow(6.4,0.25), 4.); |
| 231 |
|
| 232 |
bool TRACK__OK = false; |
| 233 |
if( |
| 234 |
// trk->chi2 >0 && |
| 235 |
trk->GetNX()>=4 && |
| 236 |
trk->GetNY()>=3 && |
| 237 |
trk->GetLeverArmX()>=5 && |
| 238 |
trk->GetRigidity()> 4. && // unit: GV |
| 239 |
trk->GetDEDX()<3. && // unit: dE/dX for Z=1 MIP |
| 240 |
trk->IsInsideCavity() && |
| 241 |
trk->chi2 < chi2m && |
| 242 |
true ) TRACK__OK = true; |
| 243 |
|
| 244 |
if( !TRACK__OK )return false; |
| 245 |
|
| 246 |
sel[2]++; |
| 247 |
|
| 248 |
|
| 249 |
//------------------------------------------------------------------ |
| 250 |
// AC |
| 251 |
//------------------------------------------------------------------ |
| 252 |
if( acl2->IsHit("CARD1-M") || acl2->IsHit("CARD1-E") )return false; |
| 253 |
if( acl2->IsHit("CARD2-M") || acl2->IsHit("CARD2-E") )return false; |
| 254 |
if( acl2->IsHit("CARD3-M") || acl2->IsHit("CARD3-E") )return false; |
| 255 |
if( acl2->IsHit("CARD4-M") || acl2->IsHit("CARD4-E") )return false; |
| 256 |
|
| 257 |
if( acl2->IsHit("CAT1-M") || acl2->IsHit("CAT1-E") )return false; |
| 258 |
if( acl2->IsHit("CAT2-M") || acl2->IsHit("CAT2-E") )return false; |
| 259 |
if( acl2->IsHit("CAT3-M") || acl2->IsHit("CAT3-E") )return false; |
| 260 |
if( acl2->IsHit("CAT4-M") || acl2->IsHit("CAT4-E") )return false; |
| 261 |
|
| 262 |
// if( acl2->IsHit("CAS1-M") || acl2->IsHit("CAS1-E") )return false; |
| 263 |
// if( acl2->IsHit("CAS2-M") || acl2->IsHit("CAS2-E") )return false; |
| 264 |
// if( acl2->IsHit("CAS3-M") || acl2->IsHit("CAS3-E") )return false; |
| 265 |
// if( acl2->IsHit("CAS4-M") || acl2->IsHit("CAS4-E") )return false; |
| 266 |
sel[3]++; |
| 267 |
|
| 268 |
|
| 269 |
|
| 270 |
//------------------------------------------------------------------ |
| 271 |
// CALORIMETER |
| 272 |
//------------------------------------------------------------------ |
| 273 |
bool NONINT=false; |
| 274 |
if( |
| 275 |
track->GetCaloTrack()->qtrack / event->GetCaloLevel2()->qtot > 0.85 && |
| 276 |
true)NONINT=true; |
| 277 |
|
| 278 |
if(!NONINT)return false; |
| 279 |
sel[4]++; |
| 280 |
|
| 281 |
|
| 282 |
return true; |
| 283 |
|
| 284 |
|
| 285 |
} |
| 286 |
//=============================================================== |
| 287 |
// Create histograms |
| 288 |
// |
| 289 |
// |
| 290 |
// |
| 291 |
// |
| 292 |
// |
| 293 |
//=============================================================== |
| 294 |
void CreateHistos( PamLevel2* event , TFile* outf ){ |
| 295 |
cout << "Creating histos..."<<endl; |
| 296 |
gROOT->cd();//create histos in memory |
| 297 |
TString hid; |
| 298 |
TString title; |
| 299 |
|
| 300 |
//-------------------------------- |
| 301 |
// histos |
| 302 |
//-------------------------------- |
| 303 |
for(Int_t i=0; i<6; i++){ |
| 304 |
hid.Form("hxytof%i",i); |
| 305 |
title.Form("Track coordinates on ToF plane %i",i); |
| 306 |
hxytof[i] = new TH2F(hid.Data(),title.Data(),6000,-30.,30.,6000,-30.,30.); |
| 307 |
} |
| 308 |
|
| 309 |
//-------------------------------- |
| 310 |
// initialization |
| 311 |
//-------------------------------- |
| 312 |
cout << "Opening DB connection..."<<endl; |
| 313 |
event->SetDBConnection(); |
| 314 |
|
| 315 |
} |
| 316 |
//=============================================================== |
| 317 |
// Retrieve histograms |
| 318 |
// |
| 319 |
// |
| 320 |
// |
| 321 |
// |
| 322 |
// |
| 323 |
//=============================================================== |
| 324 |
void RetrieveHistos(){ |
| 325 |
|
| 326 |
gROOT->cd(); |
| 327 |
TString hid; |
| 328 |
|
| 329 |
for(Int_t i=0; i<6; i++){ |
| 330 |
hid.Form("hxytof%i",i); |
| 331 |
hxytof[i] = dynamic_cast<TH2F*>(gDirectory->FindObject(hid.Data())); |
| 332 |
} |
| 333 |
|
| 334 |
} |
| 335 |
//=============================================================== |
| 336 |
// Get histograms (... non ho capito la differenza con retrieve) |
| 337 |
// |
| 338 |
// |
| 339 |
// |
| 340 |
// |
| 341 |
// |
| 342 |
//=============================================================== |
| 343 |
void GetHistos(){ |
| 344 |
|
| 345 |
TString hid; |
| 346 |
} |
| 347 |
//=============================================================== |
| 348 |
// Fill histograms (called event-by-event) |
| 349 |
// |
| 350 |
// |
| 351 |
// |
| 352 |
// |
| 353 |
// |
| 354 |
// |
| 355 |
//=============================================================== |
| 356 |
void FillHistos( PamLevel2* event ){ |
| 357 |
// ---------------------------------- |
| 358 |
// retrieve histos |
| 359 |
// ---------------------------------- |
| 360 |
RetrieveHistos(); |
| 361 |
|
| 362 |
|
| 363 |
// --------------------------------- |
| 364 |
// Retrieve the track |
| 365 |
// --------------------------------- |
| 366 |
PamTrack *track = event->GetTrack(0); |
| 367 |
TrkTrack *trk = track->GetTrkTrack(); |
| 368 |
ToFTrkVar *tof = track->GetToFTrack(); |
| 369 |
|
| 370 |
// cout << track->GetTrkTrack()->GetRigidity() << endl; |
| 371 |
// cout << track->GetToFTrack()->beta[12] << endl; |
| 372 |
// cout << track->GetCaloTrack()->qcyl << endl; |
| 373 |
|
| 374 |
// --------------------------------------------- |
| 375 |
// track coordinates esxtrapolated on tof planes |
| 376 |
// --------------------------------------------- |
| 377 |
float ytof[6],xtof[6]; |
| 378 |
for( int i=0; i<6; i++) ytof[i] = tof->ytr_tof[i]; |
| 379 |
for( int i=0; i<6; i++) xtof[i] = tof->xtr_tof[i]; |
| 380 |
|
| 381 |
for( int i=0; i<6; i++) hxytof[i]->Fill(xtof[i],ytof[i]); |
| 382 |
|
| 383 |
Dump(event->GetToFLevel2()); |
| 384 |
Dump(event->GetTrigLevel2()); |
| 385 |
UInt_t patterntrig[6]; |
| 386 |
EvaluatePatternTrig(patterntrig, event->GetToFLevel2()); |
| 387 |
cout << "TDC-Pattern: "<<hex; |
| 388 |
for(int i=0; i<6; i++)cout << setw(5) << patterntrig[i]; |
| 389 |
cout << dec << endl; |
| 390 |
cout <<"============================================="<<endl; |
| 391 |
Dump(tof); |
| 392 |
|
| 393 |
|
| 394 |
} |
| 395 |
//=============================================================== |
| 396 |
// Save histograms |
| 397 |
// |
| 398 |
// |
| 399 |
// |
| 400 |
// |
| 401 |
// |
| 402 |
//=============================================================== |
| 403 |
void SaveHistos( TFile * outf ){ |
| 404 |
|
| 405 |
TString hid; |
| 406 |
|
| 407 |
gROOT->cd(); |
| 408 |
// ---------------------------------- |
| 409 |
// retrieve histos |
| 410 |
// ---------------------------------- |
| 411 |
RetrieveHistos(); |
| 412 |
|
| 413 |
// ---------------------------------- |
| 414 |
// save on file |
| 415 |
// ---------------------------------- |
| 416 |
outf->cd(); |
| 417 |
|
| 418 |
for(Int_t i=0; i<6; i++)hxytof[i]->Write(); |
| 419 |
|
| 420 |
cout << "####------------------------------------------------------------------" <<endl; |
| 421 |
cout << "#### pre-selected :"<<tot<<endl; |
| 422 |
if(!tot)return; |
| 423 |
cout << "#### orbital :"<<sel[0]<< " (" << ((float)sel[0]/(float)tot) <<")"<<endl; |
| 424 |
cout << "#### trigger configuration :"<<sel[1]<< " (" << ((float)sel[1]/(float)tot) <<")"<<endl; |
| 425 |
cout << "#### tracker :"<<sel[2]<< " (" << ((float)sel[2]/(float)tot) <<")"<<endl; |
| 426 |
cout << "#### AC :"<<sel[3]<< " (" << ((float)sel[3]/(float)tot) <<")"<<endl; |
| 427 |
cout << "#### calo (non-int) :"<<sel[4]<< " (" << ((float)sel[4]/(float)tot) <<")"<<endl; |
| 428 |
cout << "####------------------------------------------------------------------" <<endl; |
| 429 |
|
| 430 |
} |
Parent Directory
| 
Revision Log