| 24 |
|
|
| 25 |
_whichAlg = id; |
_whichAlg = id; |
| 26 |
|
|
|
SetDebug(false); |
|
| 27 |
|
|
| 28 |
|
SetDebug(false); |
| 29 |
|
|
| 30 |
if(id == 0){ |
if(id == 0){ |
| 31 |
if(_debug){ |
if(_debug){ |
| 32 |
cout << "ExtTrkingAlg::ExtTrkingAlg("<<id<<")"<<endl; |
cout << "ExtTrkingAlg::ExtTrkingAlg("<<id<<")"<<endl; |
| 33 |
cout << "Creating array of TrkTrack objects "<<endl; |
// cout << "Creating array of TrkTrack objects "<<endl; |
| 34 |
cout << "WARNING!!! tracking not accurate!! bug not fixed yet... "<<endl; |
cout << "WARNING!!! tracking not accurate!! bug not fixed yet... "<<endl; |
| 35 |
} |
} |
| 36 |
_trkArray = new TClonesArray("TrkTrack"); |
_trkArray = new TClonesArray("TrkTrack"); |
| 55 |
}else if ( (id >= 100 && id < 144) || |
}else if ( (id >= 100 && id < 144) || |
| 56 |
(id >= 200 && id < 244) || |
(id >= 200 && id < 244) || |
| 57 |
false){ //add n calo planes |
false){ //add n calo planes |
| 58 |
if(_debug){ |
if(_debug){ |
| 59 |
cout << "ExtTrkingAlg::ExtTrkingAlg("<<id<<")"<<endl; |
cout << "ExtTrkingAlg::ExtTrkingAlg("<<id<<")"<<endl; |
| 60 |
cout << "Creating array of ExtTrack objects "<<endl; |
// cout << "Creating array of ExtTrack objects "<<endl; |
| 61 |
} |
} |
| 62 |
_trkArray = new TClonesArray("ExtTrack"); |
_trkArray = new TClonesArray("ExtTrack"); |
| 63 |
// |
// |
| 64 |
_alg_nClFixX = 3.; // algorythm parameter: n.hits required on X view |
_alg_nClFixX = 3.; // algorythm parameter: n.hits required on X view |
| 71 |
if((int)(id/200)==1)_alg_nViewCal = id - 200; // algorythm parameter: n.calorimeter planes included |
if((int)(id/200)==1)_alg_nViewCal = id - 200; // algorythm parameter: n.calorimeter planes included |
| 72 |
|
|
| 73 |
// |
// |
| 74 |
_sel_nClstrMAX = 10. ; // selection parameter: maximum number of cluster |
_sel_nClstrMAX = 15.;//10. ; // selection parameter: maximum number of cluster |
| 75 |
_sel_nPlaneXMIN = (_alg_nViewCal ? 2 : 3); // selection parameter: minimum number of hit x-views |
_sel_nPlaneXMIN = (_alg_nViewCal ? 2 : 3); // selection parameter: minimum number of hit x-views |
| 76 |
_sel_nPlaneYMIN = 2 - (int)(_alg_nViewCal+1)/2 ; // selection parameter: minimum number of hit y-views |
_sel_nPlaneYMIN = 2 - (int)(_alg_nViewCal+1)/2 ; // selection parameter: minimum number of hit y-views |
| 77 |
|
|
| 96 |
// _extTrack->SetZ(4,ZTRK5); |
// _extTrack->SetZ(4,ZTRK5); |
| 97 |
// _extTrack->SetZ(5,ZTRK6); |
// _extTrack->SetZ(5,ZTRK6); |
| 98 |
|
|
| 99 |
|
|
| 100 |
|
|
| 101 |
_caloStripRoto.clear(); |
_caloStripRoto.clear(); |
| 102 |
for(int iv=0; iv<_alg_nViewCal; iv++){ //loop over calorimater tracking planes |
for(int iv=0; iv<_alg_nViewCal; iv++){ //loop over calorimater tracking planes |
| 103 |
for(int is=0; is<9; is++){ //loop over calorimeter sensors |
for(int is=0; is<9; is++){ //loop over calorimeter sensors |
| 109 |
//------------------------------------------ |
//------------------------------------------ |
| 110 |
// read alignment parameters |
// read alignment parameters |
| 111 |
//------------------------------------------ |
//------------------------------------------ |
| 112 |
|
// printf("qui debug : %i \n",_debug); |
| 113 |
if(_debug) cout << "Reading calorimeter alignment parameters"<<endl; |
if(_debug) cout << "Reading calorimeter alignment parameters"<<endl; |
| 114 |
const char *pamca = gSystem->Getenv("PAM_CALIB"); |
const char *pamca = gSystem->Getenv("PAM_CALIB"); |
| 115 |
TString filep = "/trk-param/align_param_calo-0/CaloAlignParams.txt"; |
TString filep = "/trk-param/align_param_calo-0/CaloAlignParams.txt"; |
| 116 |
filep.Prepend(pamca); |
filep.Prepend(pamca); |
| 117 |
|
|
| 118 |
|
if(_debug)cout << " ---> "<<filep.Data()<<endl; |
| 119 |
std::fstream fs; |
std::fstream fs; |
| 120 |
fs.open (filep.Data(), std::fstream::in ); |
fs.open (filep.Data(), std::fstream::in ); |
| 121 |
do{ |
|
| 122 |
|
if(!fs.is_open()){ |
| 123 |
|
|
| 124 |
|
cout <<" WARNING!!!!!!!!!!!! error opening file "<<filep.Data()<<endl; |
| 125 |
|
cout <<" No alignment parameters loaded "<<endl; |
| 126 |
|
} |
| 127 |
|
// do{ |
| 128 |
|
while(fs.is_open() && !fs.eof()){ |
| 129 |
|
|
| 130 |
int view,plane,sensor; |
int view,plane,sensor; |
| 131 |
float par[5],err[5]; |
float par[5],err[5]; |
| 132 |
fs>>view >> plane>> sensor; |
fs>>view >> plane>> sensor; |
| 151 |
// _caloStripRoto[index].SetView(view); |
// _caloStripRoto[index].SetView(view); |
| 152 |
// _caloStripRoto[index].SetPlane(plane); |
// _caloStripRoto[index].SetPlane(plane); |
| 153 |
|
|
| 154 |
} |
}; |
| 155 |
}while(!fs.eof()); |
// }while(!fs.eof()); |
| 156 |
|
}; |
| 157 |
fs.close(); |
fs.close(); |
| 158 |
|
|
| 159 |
// ----------------------------------- |
// ----------------------------------- |
| 296 |
*/ |
*/ |
| 297 |
Bool_t ExtTrkingAlg::CheckEvent(){ |
Bool_t ExtTrkingAlg::CheckEvent(){ |
| 298 |
|
|
| 299 |
if(!_trk_l2)return true; |
if(!_trk_l1)return false; |
| 300 |
|
TrkLevel1 *trkl1 = _trk_l1;//event->GetTrkLevel1(); |
|
TrkLevel2 *trkl2 = _trk_l2;//event->GetTrkLevel2(); |
|
| 301 |
int nClstrMAX = _sel_nClstrMAX; //maximum number of cluster |
int nClstrMAX = _sel_nClstrMAX; //maximum number of cluster |
| 302 |
int nPlaneXMIN = _sel_nPlaneXMIN; |
///////////////////////////////////////////////////////////////////// |
| 303 |
int nPlaneYMIN = _sel_nPlaneYMIN; |
/// dump selection |
| 304 |
|
///////////////////////////////////////////////////////////////////// |
| 305 |
|
if(_debug){ |
| 306 |
|
cout << " n.clusters "<<trkl1->nclstr()<<endl; |
| 307 |
|
} |
| 308 |
|
if(trkl1->nclstr() > nClstrMAX) return false; |
| 309 |
///////////////////////////////////////////////////////////////////// |
///////////////////////////////////////////////////////////////////// |
| 310 |
/// count number of hit planes |
/// count number of hit planes |
| 311 |
///////////////////////////////////////////////////////////////////// |
///////////////////////////////////////////////////////////////////// |
| 312 |
|
int nPlaneXMIN = _sel_nPlaneXMIN; |
| 313 |
|
int nPlaneYMIN = _sel_nPlaneYMIN; |
| 314 |
int nHitX[] = {0,0,0,0,0,0}; |
int nHitX[] = {0,0,0,0,0,0}; |
| 315 |
int nHitY[] = {0,0,0,0,0,0}; |
int nHitY[] = {0,0,0,0,0,0}; |
| 316 |
for(int ix=0; ix<trkl2->nclsx(); ix++)nHitX[trkl2->GetSingletX(ix)->plane-1]++; |
for(int ic=0; ic<trkl1->nclstr(); ic++){ |
| 317 |
for(int iy=0; iy<trkl2->nclsy(); iy++)nHitY[trkl2->GetSingletY(iy)->plane-1]++; |
if (trkl1->GetCluster(ic)->IsX() ) nHitX[trkl1->GetCluster(ic)->GetPlane()-1]++; |
| 318 |
|
else if (trkl1->GetCluster(ic)->IsY() ) nHitY[trkl1->GetCluster(ic)->GetPlane()-1]++; |
| 319 |
|
} |
| 320 |
///////////////////////////////////////////////////////////////////// |
///////////////////////////////////////////////////////////////////// |
| 321 |
/// set minimum condition 3x+2y |
/// set minimum condition 3x+2y |
| 322 |
///////////////////////////////////////////////////////////////////// |
///////////////////////////////////////////////////////////////////// |
| 324 |
for(int ip=0; ip<6; ip++)if(nHitX[ip])nPlaneX++; |
for(int ip=0; ip<6; ip++)if(nHitX[ip])nPlaneX++; |
| 325 |
int nPlaneY=0; |
int nPlaneY=0; |
| 326 |
for(int ip=0; ip<6; ip++)if(nHitY[ip])nPlaneY++; |
for(int ip=0; ip<6; ip++)if(nHitY[ip])nPlaneY++; |
| 327 |
|
|
| 328 |
///////////////////////////////////////////////////////////////////// |
///////////////////////////////////////////////////////////////////// |
| 329 |
/// dump selection |
/// dump selection |
| 330 |
///////////////////////////////////////////////////////////////////// |
///////////////////////////////////////////////////////////////////// |
| 331 |
if(_debug){ |
if(_debug){ |
| 332 |
cout << " n.tracks "<<trkl2->GetNTracks()<<endl; |
cout << " n.x-hit planes "<<nPlaneX<<endl; |
| 333 |
cout << " n.x-singles "<<trkl2->nclsx()<<endl; |
cout << " n.y-hit planes "<<nPlaneY<<endl; |
| 334 |
cout << " n.y-singles "<<trkl2->nclsy()<<endl; |
} |
| 335 |
cout << " n.x-planes with singles "<<nPlaneX<<endl; |
if( nPlaneX<nPlaneXMIN ) return false; |
| 336 |
cout << " n.y-planes with singles "<<nPlaneY<<endl; |
if( nPlaneY<nPlaneYMIN ) return false; |
| 337 |
} |
|
| 338 |
|
|
| 339 |
|
// if(!_trk_l2)return true; |
| 340 |
|
// TrkLevel2 *trkl2 = _trk_l2;//event->GetTrkLevel2(); |
| 341 |
|
// ///////////////////////////////////////////////////////////////////// |
| 342 |
|
// /// dump selection |
| 343 |
|
// ///////////////////////////////////////////////////////////////////// |
| 344 |
|
// if(_debug){ |
| 345 |
|
// cout << " n.tracks "<<trkl2->GetNTracks()<<endl; |
| 346 |
|
// cout << " n.x-singles "<<trkl2->nclsx()<<endl; |
| 347 |
|
// cout << " n.y-singles "<<trkl2->nclsy()<<endl; |
| 348 |
|
// cout << " n.y-singles "<<trkl2->nclsy()<<endl; |
| 349 |
|
// } |
| 350 |
|
// if(trkl2->nclsy()+trkl2->nclsx() > nClstrMAX) return false; |
| 351 |
|
// ///////////////////////////////////////////////////////////////////// |
| 352 |
|
// /// count number of hit planes |
| 353 |
|
// ///////////////////////////////////////////////////////////////////// |
| 354 |
|
// int nPlaneXMIN = _sel_nPlaneXMIN; |
| 355 |
|
// int nPlaneYMIN = _sel_nPlaneYMIN; |
| 356 |
|
// int nHitX[] = {0,0,0,0,0,0}; |
| 357 |
|
// int nHitY[] = {0,0,0,0,0,0}; |
| 358 |
|
// for(int ix=0; ix<trkl2->nclsx(); ix++)nHitX[trkl2->GetSingletX(ix)->plane-1]++; |
| 359 |
|
// for(int iy=0; iy<trkl2->nclsy(); iy++)nHitY[trkl2->GetSingletY(iy)->plane-1]++; |
| 360 |
|
// ///////////////////////////////////////////////////////////////////// |
| 361 |
|
// /// set minimum condition 3x+2y |
| 362 |
|
// ///////////////////////////////////////////////////////////////////// |
| 363 |
|
// int nPlaneX=0; |
| 364 |
|
// for(int ip=0; ip<6; ip++)if(nHitX[ip])nPlaneX++; |
| 365 |
|
// int nPlaneY=0; |
| 366 |
|
// for(int ip=0; ip<6; ip++)if(nHitY[ip])nPlaneY++; |
| 367 |
|
// ///////////////////////////////////////////////////////////////////// |
| 368 |
|
// /// dump selection |
| 369 |
|
// ///////////////////////////////////////////////////////////////////// |
| 370 |
|
// if(_debug){ |
| 371 |
|
// cout << " n.x-planes with singles "<<nPlaneX<<endl; |
| 372 |
|
// cout << " n.y-planes with singles "<<nPlaneY<<endl; |
| 373 |
|
// } |
| 374 |
//////////////////////////////////////////////////////////////////////// |
//////////////////////////////////////////////////////////////////////// |
| 375 |
// |
// |
| 376 |
// check condition for retracking |
// check condition for retracking |
| 377 |
// |
// |
| 378 |
//////////////////////////////////////////////////////////////////////// |
//////////////////////////////////////////////////////////////////////// |
| 379 |
if(trkl2->GetNTracks()==0 && nPlaneX<nPlaneXMIN) return false; |
// if(trkl2->GetNTracks()==0 && nPlaneX<nPlaneXMIN) return false; |
| 380 |
if(trkl2->GetNTracks()==0 && nPlaneY<nPlaneYMIN) return false; |
// if(trkl2->GetNTracks()==0 && nPlaneY<nPlaneYMIN) return false; |
|
if(trkl2->nclsy()+trkl2->nclsx() > nClstrMAX) return false; |
|
| 381 |
|
|
| 382 |
return true; |
return true; |
| 383 |
|
|
| 1062 |
Float_t mip = l1->DecodeEstrip(ih,view,plane,strip); |
Float_t mip = l1->DecodeEstrip(ih,view,plane,strip); |
| 1063 |
|
|
| 1064 |
// if(strip == 0 || strip == 95)cout <<" strip "<<strip<<" OK "<<endl; |
// if(strip == 0 || strip == 95)cout <<" strip "<<strip<<" OK "<<endl; |
| 1065 |
if(strip < 0 || strip > 95)cout <<" strip "<<strip<<" AHI AHI "<<endl; |
if((strip < 0 || strip > 95) && _debug)cout <<" strip "<<strip<<" AHI AHI "<<endl; |
| 1066 |
|
|
| 1067 |
if(view<0)continue; |
if(view<0)continue; |
| 1068 |
if(plane<0)continue; |
if(plane<0)continue; |
| 1419 |
// cout << " extst.Set("<<(int)isY<<", "<<plane<<", "<<cluster.coordPU<<", "<<sensor<<") "<<endl; |
// cout << " extst.Set("<<(int)isY<<", "<<plane<<", "<<cluster.coordPU<<", "<<sensor<<") "<<endl; |
| 1420 |
// extst.Set((int)isY, plane, cluster.coordPU, sensor); |
// extst.Set((int)isY, plane, cluster.coordPU, sensor); |
| 1421 |
|
|
| 1422 |
if( sisensor != extst.GetSiSensor() )cout << " extst.GetSiSensor() "<<extst.GetSiSensor()<<" != "<<sisensor<<endl; |
if( sisensor != extst.GetSiSensor() && _debug )cout << " extst.GetSiSensor() "<<extst.GetSiSensor()<<" != "<<sisensor<<endl; |
| 1423 |
if( plane != extst.GetPlane() )cout << " extst.GetPlane() "<<extst.GetPlane()<<" != "<<plane<<endl; |
if( plane != extst.GetPlane() && _debug )cout << " extst.GetPlane() "<<extst.GetPlane()<<" != "<<plane<<endl; |
| 1424 |
if( (int)isY != extst.GetView() )cout << " extst.GetView() "<<extst.GetView()<<" != "<<isY<<endl; |
if( (int)isY != extst.GetView() && _debug )cout << " extst.GetView() "<<extst.GetView()<<" != "<<isY<<endl; |
| 1425 |
// cout << " extst.SetStrip("<<cluster.coordPU<<") "<<endl; |
// cout << " extst.SetStrip("<<cluster.coordPU<<") "<<endl; |
| 1426 |
extst.SetStrip(cluster.coordPU); |
extst.SetStrip(cluster.coordPU); |
| 1427 |
|
|
| 1525 |
// if(isY) VMS[nW] = rms*rms; |
// if(isY) VMS[nW] = rms*rms; |
| 1526 |
// if(_debug)cout <<endl<<view<<" nW "<<nW<<" res(MS) "<<rms<<" factor "<<factor; |
// if(_debug)cout <<endl<<view<<" nW "<<nW<<" res(MS) "<<rms<<" factor "<<factor; |
| 1527 |
|
|
| 1528 |
if(view>=3)cout << " --> NB error computation not accurate "; |
if(view>=3 && _debug)cout << " --> NB error computation not accurate "; |
| 1529 |
res_ms = rms; |
res_ms = rms; |
| 1530 |
|
|
| 1531 |
} |
} |
| 2929 |
int ifail=0; |
int ifail=0; |
| 2930 |
|
|
| 2931 |
|
|
| 2932 |
|
// trackCand.Dump(); |
| 2933 |
|
|
| 2934 |
trackCand.ResetFit(); |
trackCand.ResetFit(); |
| 2935 |
trackCand.Fit(0.,ifail,0); |
trackCand.Fit(0.,ifail,0); |
| 2936 |
if(ifail!=0)trackCand.ResetFit(); |
if(ifail!=0)trackCand.ResetFit(); |
| 2937 |
|
|
| 2938 |
// trackCand.Dump(); |
// trackCand.Dump(); |
| 2939 |
|
|
| 2940 |
|
|
| 2941 |
// ----------------------------------- |
// ----------------------------------- |
| 2977 |
if(_debug)cout << "n.candidates "<<trackCandidates[mapIndex].size()<<endl; |
if(_debug)cout << "n.candidates "<<trackCandidates[mapIndex].size()<<endl; |
| 2978 |
|
|
| 2979 |
|
|
|
|
|
| 2980 |
//================================================================= |
//================================================================= |
| 2981 |
// ---------------------------------------------------------------- |
// ---------------------------------------------------------------- |
| 2982 |
// CALORIMETER |
// CALORIMETER |