MIPP(E907): SPFit Class Reference

SPFit::SPFit ( const char * version )

Definition at line 35 of file SPFit.cxx.

References JobCModule::CheckInit(), fVtxResidNt, and JobCModule::SetCfgVersion().

00035                                 :
00036   JobCModule("SPFit"),
00037   fResidNt(0),
00038   fTrackNt(0),
00039   fVtxResidNt(0),
00040   fTPCResCorMap(0),
00041   fTPCResCorMaker(0),
00042   fGeo(0)
00043 {
00044   this->SetCfgVersion(version);
00045   CheckInit();
00046 
00047   fVtxResidNt = new TNtuple("fVtxResidNt","SPFit Vertexing Ntuple",
00048                 "bx:by:tx:ty:ztgt,mom");  
00049 
00050 }

SPFit::~SPFit ( )

Definition at line 54 of file SPFit.cxx.

References fResidNt, fTPCResCorMaker, fTrackNt, and fVtxResidNt.

00055 {
00056   if (fResidNt) { fResidNt->Write(); delete fResidNt; }
00057   if (fTrackNt) { fTrackNt->Write(); delete fTrackNt; }
00058   if (fVtxResidNt) { fVtxResidNt->Write(); delete fVtxResidNt; }
00059   if (fTPCResCorMaker) {
00060     delete fTPCResCorMaker; fTPCResCorMaker = 0;
00061   }
00062 }

void SPFit::EndRun	(	int	run,
		int	subrun
	)			`[virtual]`

Reimplemented from JobCModule.

Definition at line 116 of file SPFit.cxx.

References TPCResCorMaker::Finalize(), fTPCResCorMaker, fTPCResCorMap, fUpdateTPCResCor, TPCResCorMaker::WriteHisto(), and TPCResCorMap::WriteTableToFile().

00117 {
00118   if (fUpdateTPCResCor) {
00119     fTPCResCorMaker->Finalize();
00120     fTPCResCorMap->WriteTableToFile("spfit-rescor.txt");
00121     fTPCResCorMaker->WriteHisto("tpc-res-cor.root");
00122   }
00123 }

void SPFit::NewRun	(	int	run,
		int	subrun
	)			`[virtual]`

Reimplemented from JobCModule.

Definition at line 66 of file SPFit.cxx.

References RIBfield::AvgJGG(), RIBfield::AvgRos(), fGeo, fJGGField, fROSYField, fTPCResCorMaker, fTPCResCorMap, fTPCResCorMapSrc, fUpdateTPCResCor, fUseTPCResCor, GMIPPGeo::Instance(), and RIBfield::Instance().

00067 {
00068   RIBfield& rib = RIBfield::Instance(run);
00069 
00070   fGeo = & GMIPPGeo::Instance(run);
00071 
00072   fJGGField  = rib.AvgJGG(2);
00073   fROSYField = rib.AvgRos(2);
00074 
00075   // Trap Nan's
00076   if (fJGGField!=fJGGField)   {
00077     fJGGField  =  0.70;
00078     std::cout << __FILE__ << ":" << __LINE__ 
00079           << " Failed to load JGG field for run " << run
00080           << " Using a default value of " << fJGGField << std::endl;
00081   }
00082   if (fROSYField!=fROSYField) {
00083     fROSYField = -0.54;
00084     std::cout << __FILE__ << ":" << __LINE__ 
00085           << " Failed to load ROSY field for run " << run
00086           << " Using a default value of " << fROSYField << std::endl;
00087   }
00088 
00089   // If we've been asked to use the TPC residual corrections, build
00090   // the map
00091   if (fUseTPCResCor) {
00092     if (fTPCResCorMap) { 
00093       delete fTPCResCorMap;
00094       fTPCResCorMap = 0; 
00095     }
00096     fTPCResCorMap = new TPCResCorMap(run,fTPCResCorMapSrc.c_str());
00097     std::cerr << "Read file " << fTPCResCorMapSrc << std::endl;
00098     std::cerr << " " << fTPCResCorMap->Nx() 
00099           << " " << fTPCResCorMap->Ny()
00100           << " " << fTPCResCorMap->Npadrow() << std::endl;
00101     
00102     // If we've been asked to update the residual corrections, build
00103     // the class that computes the updates
00104     if (fUpdateTPCResCor) {
00105       if (fTPCResCorMaker) {
00106     delete fTPCResCorMaker;
00107     fTPCResCorMaker = 0; 
00108       }
00109       fTPCResCorMaker = new TPCResCorMaker(fTPCResCorMap);
00110     }
00111   }
00112 }

JobCResult SPFit::Reco ( EDMEventHandle & evt ) [virtual]

Reimplemented from JobCModule.

Definition at line 127 of file SPFit.cxx.

References SPTrackNt::evt, SPResidNt::evt, fDebug, fDoSimpleVtx, fFillResidNt, fFillTrackNt, fGeo, fJGGField, fPmax, fPmin, fResidNt, fROSYField, fTPCResCorMaker, fTPCResCorMap, fTPCXWeight, fTPCYWeight, fTrackNt, fVtxResidNt, EDMDataBucket::Get(), EDMEventHandle::Header(), SPTrackNt::idc, SPResidNt::idc, SPTrackNt::isci, SPResidNt::isci, JobCModule::kFailed, JobCModule::kPassed, GMIPPGeo::kThin, EDMDataBucket::kVariableSizeObject, mom, SPTrackNt::ntrk, EDMDataBucket::Put(), EDMEventHandle::Reco(), SPTrackNt::run, SPResidNt::run, EDMEventHeader::Run(), GMIPPGeo::Target(), GMIPPGeo::TargetType(), SPTrackNt::trk, SPResidNt::trk, x, y, GTargetGeo::ZMax(), and GTargetGeo::ZMin().

00128 {
00129   unsigned int i, j;
00130 
00131   // Grab the track clusters out of the event
00132   std::vector<const SPTrkClust*> trkseg;
00133   try                    { evt.Reco().Get("./spfit",trkseg); }
00134   catch (EDMException e) { 
00135     if (fDebug>0) {
00136       std::cerr << "SPFit: Failed to find track segments." << std::endl;
00137     }
00138     return JobCModule::kFailed; 
00139   }
00140 
00141   int irun = evt.Header().Run();
00142 
00143   int maxTPCid = 0;
00144   // Find largest TPC track ID number:
00145   for (i=0; i<trkseg.size(); ++i) {
00146     if (trkseg[i]->TPC()) {
00147       int id = trkseg[i]->TPC()->Id();
00148       if (id > maxTPCid) maxTPCid = id;
00149     }
00150   }
00151   
00152   // Fit the track segments associated with the track during previous
00153   // stage
00154   std::vector<SPTrack*> sptrk;
00155   for (i=0; i<trkseg.size(); ++i) {
00156     SPTrack* t = new SPTrack(*trkseg[i]);
00157     
00158     t->SetWeights(fTPCXWeight,fTPCYWeight);
00159     t->SetBFields(fJGGField,fROSYField);
00160     t->SeedFit(fGeo);
00161     if (fTPCResCorMap) t->ApplyTPCResCor(fTPCResCorMap);
00162 
00163     // If the momentum is within range, try to improve the
00164     // fit. Otherwise, just use the seed results.
00165     if (t->GetMomentum()>=fPmin && t->GetMomentum()<=fPmax) t->Fit(fGeo);
00166     else                                                    t->Fit(fGeo,0);
00167 
00168     if (trkseg[i]->TPC()) t->SetId(trkseg[i]->TPC()->Id());    
00169     else                  t->SetId(++maxTPCid);
00170     sptrk.push_back(t);
00171   }
00172 
00173   // Preen the list of tracks to remove duplicate tracks
00174   if (sptrk.size()>1) {
00175     for (i=0; i<sptrk.size()-1; ++i) {
00176       for (j=i+1; j<sptrk.size(); ++j) {
00177     if (sptrk[i] == 0)              continue;
00178     if (sptrk[j] == 0)              continue;
00179     if (sptrk[i]->Clust()->TPC()==0) continue;
00180     if (sptrk[j]->Clust()->TPC()==0) continue;
00181     if (sptrk[i]->Clust()->TPC()->Id() == 
00182         sptrk[j]->Clust()->TPC()->Id()) {
00183       // Tracks have same TPC tracks associated -- keep the one
00184       // that fits better
00185       if (sptrk[j]->GoF() > sptrk[i]->GoF()) {
00186         if (fDebug>0) {
00187           std::cout << "SPFit: Track " 
00188             << j << "/" << sptrk[j]->GoF() 
00189             << " beats "
00190             << i << "/" << sptrk[i]->GoF()
00191             << std::endl;
00192         }
00193         delete sptrk[i]; sptrk[i] = 0;
00194       }
00195       else {
00196         if (fDebug>0) {
00197           std::cout << "SPFit: Track " 
00198             << i << "/" << sptrk[i]->GoF()
00199             << " beats "
00200             << j << "/" << sptrk[j]->GoF()
00201             << std::endl;
00202         }
00203         delete sptrk[j]; sptrk[j] = 0;
00204       }
00205     } // Tracks share same TPC segment
00206       } // Loop on track j
00207     } // Loop on track i
00208   } // Track list has non-zero size
00209   
00210   vector<RBKTrack*> rbtlist(0);
00211   for (i=0; i<sptrk.size(); ++i) {
00212     SPTrack* t = sptrk[i];
00213     RBKTrack* rbt = 0;
00214     if (t) rbt = t->AsRBKTrack();
00215     rbtlist.push_back(rbt);
00216   }
00217   
00218   // Now do simple vertexing...
00219 
00220   if (fDoSimpleVtx && rbtlist.size() > 1) {
00221     double bx[3], dbx[3], mom[3], ztgt;
00222     double x, y;
00223     
00224     vector<const TrackSeg*> bctrk;
00225     try { 
00226       evt.Reco().Get("./TrkRBaseBC",bctrk); 
00227     }
00228     catch (EDMException e) { }
00229     
00230     if (bctrk.size()==1) {
00231       
00232       if (fGeo->TargetType() == GMIPPGeo::kThin) {
00233     cout << "SPFit: trying to find vertex for thin target!" << endl;
00234     ztgt = (fGeo->Target().ZMax() + fGeo->Target().ZMin())/2.;
00235     ztgt += 2.;
00236     cout << "ztgt = " << ztgt << endl;
00237 
00238     bctrk[0]->GetXY(ztgt,bx,dbx);
00239     
00240     cout << "bx = " << bx[0] << endl;
00241     cout << "by = " << bx[1] << endl;
00242 
00243     //    vector<double>dxtgt(0);
00244     //    vector<double>dytgt(0);
00245     //    vector<RBKTrack*>trklist(0);
00246     
00247     for (i=0; i<rbtlist.size(); ++i) {
00248       RBKTrack* rbt = rbtlist[i];
00249       if (! rbt) continue;
00250       if (! rbt->GetXY(ztgt, &x, &y)) continue;
00251 
00252       cout << "x = " << x << endl;
00253       cout << "y = " << y << endl;
00254       rbt->GetPatZ(ztgt,mom);
00255       double ptot = sqrt(mom[0]*mom[0] + mom[1]*mom[1] + mom[2]*mom[2]);
00256 
00257       fVtxResidNt->Fill(bx[0],bx[1],x,y,ztgt,ptot);
00258       
00259       //        dxtgt.push_back(x-bx);
00260       //        dytgt.push_back(y-by);  
00261       //        trklist.push_back(rbt);
00262     }
00263       }
00264     }
00265   }
00266   
00267   // Store tracks and fill ntuple summaries
00268   for (i=0; i<rbtlist.size(); ++i) {
00269 
00270     RBKTrack* rbt = rbtlist[i];
00271     if (!rbt) continue;
00272     if (!sptrk[i]) continue;
00273 
00274     rbt->SetId(sptrk[i]->GetId());
00275     
00276     evt.Reco().Put(*rbt, "./spfit",
00277            EDMDataBucket::kVariableSizeObject);
00278 
00279     if (fDebug>0) sptrk[i]->PrintSummary();
00280     
00281     // Fill track and residual ntuples if requested
00282     if (fFillResidNt>0) {
00283       static struct SPResidNt* sResidNt = 0;
00284       if (sResidNt==0) sResidNt = new SPResidNt;
00285       if (fResidNt==0) {
00286     fResidNt = new TTree("fSPResidNt","SPFit Residuals");
00287     fResidNt->Branch("fSPResidNt","SPResidNt",&sResidNt);
00288       }
00289       sResidNt->run  = irun;
00290       sResidNt->evt  = evt.Header().Event();
00291       sResidNt->trk  = i;
00292       sResidNt->idc  = 0;
00293       sResidNt->isci = 0;
00294       sptrk[i]->FillResidualNt(*fResidNt,*sResidNt);
00295     }
00296 
00297     // Fill the TPC residual tables, if requested
00298     if (fTPCResCorMaker) sptrk[i]->FillTPCResCor(fTPCResCorMaker);
00299 
00300     if (fFillTrackNt>0) {
00301       static struct SPTrackNt* sTrackNt = 0;
00302       if (sTrackNt==0) sTrackNt = new SPTrackNt;
00303       if (fTrackNt==0) {
00304     fTrackNt = new TTree("fSPTrackNt","SPFit Tracks");
00305     fTrackNt->Branch("fSPTrackNt","SPTrackNt",&sTrackNt);
00306       }
00307       sTrackNt->run = irun;
00308       sTrackNt->evt = evt.Header().Event();
00309       sTrackNt->trk = i;
00310       sTrackNt->idc = 0;
00311       sTrackNt->isci = 0;
00312       sTrackNt->ntrk = trkseg.size();
00313       sptrk[i]->FillTrackNt(*sTrackNt);
00314       fTrackNt->Fill();
00315     }
00316   }
00317 
00318   // Clean up
00319   for (i=0; i<sptrk.size(); ++i) {
00320     if (sptrk[i]!=0) { delete sptrk[i]; sptrk[i] = 0; }
00321   }
00322   
00323   return JobCModule::kPassed;
00324 }

void SPFit::Update ( const CfgConfig & c ) [virtual]

Implements CfgObserver.

Definition at line 328 of file SPFit.cxx.

References fDebug, fDoSimpleVtx, fFillResidNt, fFillTrackNt, JobCModule::fIsInit, fPmax, fPmin, fTPCResCorMapSrc, fTPCXWeight, fTPCYWeight, fUpdateTPCResCor, and fUseTPCResCor.

00329 {
00330   c("Debug").        Get(fDebug);
00331   c("FillTrackNt").  Get(fFillTrackNt);
00332   c("FillResidNt").  Get(fFillResidNt);
00333   c("TPCXWeight").   Get(fTPCXWeight);
00334   c("TPCYWeight").   Get(fTPCYWeight);
00335   c("Pmin").         Get(fPmin);
00336   c("Pmax").         Get(fPmax);
00337   c("DoSimpleVtx").  Get(fDoSimpleVtx);
00338 
00339   c("UseTPCResCor").   Get(fUseTPCResCor);
00340   c("UpdateTPCResCor").Get(fUpdateTPCResCor);
00341   c("TPCResCorMapSrc").Get(fTPCResCorMapSrc);
00342 
00343   fIsInit = true;
00344 }

int SPFit::fDebug [private]

Debug level 0=none, >0 gets more verbose.

Definition at line 40 of file SPFit.h.

Referenced by Reco(), and Update().

int SPFit::fDoSimpleVtx [private]

flag to do simple vertexing

Definition at line 43 of file SPFit.h.

Referenced by Reco(), and Update().

int SPFit::fFillResidNt [private]

1 = yes, 0 = no

Definition at line 42 of file SPFit.h.

Referenced by Reco(), and Update().

int SPFit::fFillTrackNt [private]

1 = yes, 0 = no

Definition at line 41 of file SPFit.h.

Referenced by Reco(), and Update().

GMIPPGeo* SPFit::fGeo [private]

Definition at line 57 of file SPFit.h.

Referenced by NewRun(), and Reco().

double SPFit::fJGGField [private]

Average field in JGG for run.

Definition at line 52 of file SPFit.h.

Referenced by NewRun(), and Reco().

float SPFit::fPmax [private]

Upper limit of momentum of tracks to fit.

Definition at line 55 of file SPFit.h.

Referenced by Reco(), and Update().

float SPFit::fPmin [private]

Lower limit of momentum of tracks to fit.

Definition at line 54 of file SPFit.h.

Referenced by Reco(), and Update().

TTree* SPFit::fResidNt [private]

Ntuple of residual information.

Definition at line 33 of file SPFit.h.

Referenced by Reco(), and ~SPFit().

double SPFit::fROSYField [private]

Average field in ROSY for run.

Definition at line 53 of file SPFit.h.

Referenced by NewRun(), and Reco().

TPCResCorMaker* SPFit::fTPCResCorMaker [private]

Definition at line 37 of file SPFit.h.

Referenced by EndRun(), NewRun(), Reco(), and ~SPFit().

TPCResCorMap* SPFit::fTPCResCorMap [private]

Definition at line 36 of file SPFit.h.

Referenced by EndRun(), NewRun(), and Reco().

std::string SPFit::fTPCResCorMapSrc [private]

Source of residual correction map.

Definition at line 49 of file SPFit.h.

Referenced by NewRun(), and Update().

float SPFit::fTPCXWeight [private]

Scale on TPC x measurements.

Definition at line 44 of file SPFit.h.

Referenced by Reco(), and Update().

float SPFit::fTPCYWeight [private]

Scale on TPC y measurements.

Definition at line 45 of file SPFit.h.

Referenced by Reco(), and Update().

TTree* SPFit::fTrackNt [private]

Ntuple for tracks information.

Definition at line 34 of file SPFit.h.

Referenced by Reco(), and ~SPFit().

int SPFit::fUpdateTPCResCor [private]

Should we update the correction table?

Definition at line 48 of file SPFit.h.

Referenced by EndRun(), NewRun(), and Update().

int SPFit::fUseTPCResCor [private]

Apply residual corrections to TPC hits?

Definition at line 47 of file SPFit.h.

Referenced by NewRun(), and Update().

TNtuple* SPFit::fVtxResidNt [private]

Simple ntuple for vertex studies.

Definition at line 35 of file SPFit.h.

Referenced by Reco(), SPFit(), and ~SPFit().

SPFit Class Reference

Public Member Functions

Private Attributes

Detailed Description

Constructor & Destructor Documentation

Member Function Documentation

Member Data Documentation


Public Member Functions
	SPFit (const char *version)
	~SPFit ()
void	NewRun (int run, int subrun)
void	EndRun (int run, int subrun)
JobCResult	Reco (EDMEventHandle &evt)
void	Update (const CfgConfig &c)
Private Attributes
TTree *	fResidNt
	Ntuple of residual information.
TTree *	fTrackNt
	Ntuple for tracks information.
TNtuple *	fVtxResidNt
	Simple ntuple for vertex studies.
TPCResCorMap *	fTPCResCorMap
TPCResCorMaker *	fTPCResCorMaker
int	fDebug
	Debug level 0=none, >0 gets more verbose.
int	fFillTrackNt
	1 = yes, 0 = no
int	fFillResidNt
	1 = yes, 0 = no
int	fDoSimpleVtx
	flag to do simple vertexing
float	fTPCXWeight
	Scale on TPC x measurements.
float	fTPCYWeight
	Scale on TPC y measurements.
int	fUseTPCResCor
	Apply residual corrections to TPC hits?
int	fUpdateTPCResCor
	Should we update the correction table?
std::string	fTPCResCorMapSrc
	Source of residual correction map.
double	fJGGField
	Average field in JGG for run.
double	fROSYField
	Average field in ROSY for run.
float	fPmin
	Lower limit of momentum of tracks to fit.
float	fPmax
	Upper limit of momentum of tracks to fit.
GMIPPGeo *	fGeo