|
//____________________________________________________________________ // // // //____________________________________________________________________ // // $Id: BrTpcTrackFollowModule.cxx,v 1.5 2002/01/03 19:53:42 cholm Exp $ // $Author: cholm $ // $Date: 2002/01/03 19:53:42 $ // $Copyright: (C) 2001 BRAHMS Collaboration <brahmlib@rhic.bnl.gov> // #ifndef BRAT_BrTpcTrackFollowModule #include "BrTpcTrackFollowModule.h" #endif #ifndef ROOT_TDirectory #include "TDirectory.h" #endif #ifndef BRAT_BrGeometryDbManager #include "BrGeometryDbManager.h" #endif #ifndef BRAT_BrTableNames #include "BrTableNames.h" #endif #ifndef WIN32 #include <iostream> #include <iomanip> #else #include <iostream.h> #include <iomanip.h> #endif //____________________________________________________________________ ClassImp(BrTpcTrackFollowModule); //____________________________________________________________________ BrTpcTrackFollowModule::BrTpcTrackFollowModule() : BrTpcTrackModule() { // Default constructor. DO NOT USE SetState(kSetup); SetMaxRowsMissed(); SetSearchWidthX(); SetSearchWidthY(); } //____________________________________________________________________ BrTpcTrackFollowModule::BrTpcTrackFollowModule(const Char_t* name, const Char_t* title) : BrTpcTrackModule(name, title) { // Named Constructor SetState(kSetup); SetMaxRowsMissed(); SetSearchWidthX(); SetSearchWidthY(); } //____________________________________________________________________ void BrTpcTrackFollowModule::DefineHistograms() { // Define histograms. They are: // <fill in here> if (GetState() != kInit) { Stop("DefineHistograms", "Must be called after Init"); return; } TDirectory* saveDir = gDirectory; Char_t dirName[32]; sprintf(dirName,"TrackFollow_%s",GetName()); TDirectory* histDir = gDirectory->mkdir(dirName); if (!histDir) { Warning("DefineHistograms","Could not create histogram subdirectory"); return; } histDir->cd(); Char_t HistName[48]; Char_t HistTitle[48]; const Int_t nRows = fParams->GetNumberOfRows(); sprintf(HistName,"hHitSearchDx_%s",GetName()); sprintf(HistTitle,"Hit - Proj search deviation (X) %s",GetName()); hHitSearchDx = new TH1F(HistName, HistTitle, 100, -2.0,2.0); sprintf(HistName,"hHitRowSearchDx_%s",GetName()); sprintf(HistTitle,"Hit - Proj search deviation vs Row (X) %s",GetName()); hHitRowSearchDx = new TH2F(HistName, HistTitle, 100, -2.0, 2.0, nRows, 1, nRows+1); sprintf(HistName,"hHitSearchDy_%s",GetName()); sprintf(HistTitle,"Hit - Proj search deviation (Y) %s",GetName()); hHitSearchDy = new TH1F(HistName, HistTitle, 100, -2.0,2.0); sprintf(HistName,"hHitRowSearchDy_%s",GetName()); sprintf(HistTitle,"Hit - Proj search deviation vs Row (Y) %s",GetName()); hHitRowSearchDy = new TH2F(HistName, HistTitle, 100, -2.0, 2.0, nRows, 1, nRows+1); sprintf(HistName,"hHitAcceptedDx_%s",GetName()); sprintf(HistTitle,"Hit - Proj deviation (X) %s",GetName()); hHitAcceptedDx = new TH1F(HistName, HistTitle, 100, -2.0,2.0); sprintf(HistName,"hHitRowAcceptedDx_%s",GetName()); sprintf(HistTitle,"Hit - Proj deviation vs Row (X) %s",GetName()); hHitRowAcceptedDx = new TH2F(HistName, HistTitle, 100, -2.0, 2.0, nRows, 1, nRows+1); sprintf(HistName,"hHitAcceptedDy_%s",GetName()); sprintf(HistTitle,"Hit - Proj deviation (Y) %s",GetName()); hHitAcceptedDy = new TH1F(HistName, HistTitle, 100, -2.0,2.0); sprintf(HistName,"hHitRowAcceptedDy_%s",GetName()); sprintf(HistTitle,"Hit - Proj deviation vs Row (Y) %s",GetName()); hHitRowAcceptedDy = new TH2F(HistName, HistTitle, 100, -2.0, 2.0, nRows, 1, nRows+1); // also needs row based histograms of search BrTpcTrackModule::DefineHistograms(); gDirectory = saveDir; } //____________________________________________________________________ void BrTpcTrackFollowModule::Init() { // Job-level initialisation // Init BrTpcTrackModule as well SetState(kInit); BrTpcTrackModule::Init(); BrGeometryDbManager* geomDb = BrGeometryDbManager::Instance(); fTpcVolume = (BrDetectorVolume*)geomDb-> GetDetectorVolume("BrDetectorVolume", (char*)GetName()); } //____________________________________________________________________ void BrTpcTrackFollowModule::Event(BrEventNode* inNode, BrEventNode* outNode) { // Event method. // First local tables and clonesarrays are cleared. // Then hits are filled in the clonesarray for each row.. // The algorithm loops over hits in the last row. looking at the // next active row in a search window that covers the front of the // TPC the following windows is determined by the fSearchWidthX - and // Y - which are angles. // // After the tracking the outputnode is filled with the best tracks // // If the tracking fails because of too many local track candidates // (> 5000) then the Status is set to 1 and can be accessed with // GetStatus. 0 is okay. SetState(kEvent); // get the input node fHitTable = inNode->GetDataTable(Form("%s %s", BRTABLENAMES kTpcHit, GetName())); if(fHitTable == 0) { if(DebugLevel() > 0 ) Warning("Event", "BrDataTable %s %s not found in this event!", BRTABLENAMES kTpcHit, GetName()); return; } // Set status to ok fStatus = 0; if(DebugLevel()) cout << endl << Form("Entering BrTpcTrackFollowModule %s Event", GetName()) << endl; BrTpcTrackModule::Clear(); BrTpcTrackModule::SortHitsIntoRowTable(); FindTracks(); if( HistOn()) hStatus->Fill(fStatus); if (fStatus) { return; } BrTpcTrackModule::FillTrackGroups(); BrTpcTrackModule::SelectBestTracks(outNode); } //_________________________________________________________________________ void BrTpcTrackFollowModule::FindTracks(){ // // Loop over rows to find tracks. // Going from the back of the Tpc to the front if(DebugLevel()) { cout << "Entering FindTracks" << endl; Print(); } const Int_t nrow = fParams->GetNumberOfActiveRows(); const Int_t lastrow = 1; for(Int_t i = 0; i <= fMaxRowsMissed; i++) { // cout << endl << "Now doing active row" << nrow-i // << "/ " << nrow << endl; FindTracksFromRow(nrow-i, lastrow); if (fStatus) return; } } //_________________________________________________________________________ void BrTpcTrackFollowModule::FindTracksFromRow(const Int_t istart, const Int_t ilast) { // This is the method that loops over the hits in a single row // It calls BrTrackSingleHit for each hit if(DebugLevel()>1){ cout << "Entering FindLocalTracksFromRow" << endl; } const Int_t rowNo = fParams->GetActiveRowNumber(istart); const Int_t numHits = fHitsInRow[rowNo-1].GetEntriesFast(); const Int_t nTracksBefore = GetNTrackCandidates(); if(DebugLevel()>2){ cout << "Hits for row no. " << rowNo << " : " << numHits << endl; } for(Int_t i = 0; i < numHits; i++) { // cout << "Doing tracking for hit : " << i << endl; BrTpcHit *hit = (BrTpcHit*)fHitsInRow[rowNo-1].At(i); // hit->Print(); cout << endl; // Check that hit have not already been used if(hit->GetUsed() > 0) continue; BrTpcTrackCandidate* track = BrTpcTrackModule::AddNewTrackCandidate(); // Test that limit is not crossed if(fStatus) return; BrVector3D pos(hit->GetPos()); track->SetPos(pos); track->SetLastRow(rowNo); track->SetStatus(BrTrackCandidate::kOk); track->AddHit(hit); if(DebugLevel()>2){ cout << "Adding new track (No. " << BrTpcTrackModule::GetNTrackCandidates() <<") for : " << endl; hit->Print(); } // do tracking with that track candidate TrackSingleHit(istart, ilast); if(fStatus) return; } // After we have done tracking for the whole row mark new tracks as used const Int_t nTracksAfter = GetNTrackCandidates(); for(Int_t j = nTracksBefore; j < nTracksAfter; j++ ) { BrTpcTrackCandidate *track = BrTpcTrackModule::GetTrackCandidateAt(j); track->IncrementHitsUsed(); } } //_________________________________________________________________________ void BrTpcTrackFollowModule::TrackSingleHit(const Int_t istart, const Int_t ilast) { // This routine tracks a single hit through all the remaining rows // and forks the track candidate into more if more than 1 hit is // found in the search window // Number of tracks at the start of this routine const Int_t numTracksBefore = GetNTrackCandidates(); // used later TObjArray FoundHits; // since we look ahead we should stop before the last active row // (= active row 2) for(Int_t irow = istart; irow >= ilast+1; irow--) { const Int_t numTracksCurrent = GetNTrackCandidates(); const Int_t nextRow = fParams->GetActiveRowNumber(irow-1); const Float_t nextRowZ = fParams->GetRowPosition(nextRow); if(DebugLevel() > 2) { cout << "next row " << setw(2) << nextRow << "," << " Number of candidates : " << numTracksCurrent << endl; } for(Int_t i = numTracksBefore-1; i < numTracksCurrent; i++) { BrTpcTrackCandidate* track = (BrTpcTrackCandidate*)GetTrackCandidateAt(i); // If bad continue if(track->IsBad()) continue; // Difference between row from last hit and nextRow // Used for scaling search windows const Float_t dRow = track->GetLastRow() - nextRow; if(dRow <= 0) { Warning("TrackSingleHit", "This should never happen!!!!!!"); } // project to nextRow and effective search windows Float_t proj[2]; Float_t dx1, dx2, dy1, dy2; if(track->GetNhit() == 1) { // best guess is right ahead proj[0] = track->GetPos()[0]; proj[1] = track->GetPos()[1]; const Float_t *const size = fTpcVolume->GetSize(); const Float_t dz = (track->GetPos()[2] + size[ 2 ]/Float_t(2)) / fParams->GetRowDistance(); dx1 = track->GetPos()[0] + size[0]/Float_t(2); dx2 = size[ 0 ] - dx1; // i.e., the rest dx1 = TMath::Max(dx1, fSearchWidthX); dx2 = TMath::Max(dx2, fSearchWidthX); dx1 *= dRow / dz; dx2 *= dRow / dz; dy1 = track->GetPos()[1] + size[1]/Float_t(2); dy2 = size[ 1 ] - dy1; // i.e., the rest dy1 = TMath::Max(dy1, fSearchWidthY); dy2 = TMath::Max(dy2, fSearchWidthY); dy1 *= dRow / dz; dy2 *= dRow / dz; } else { // projection track->GetXYPositionAtZ(proj[0], proj[1], nextRowZ); // use symmetric search windows dx1 = fSearchWidthX * dRow; dx2 = dx1; dy1 = fSearchWidthY * dRow; dy2 = dy1; } FindHitsNextRow(nextRow, proj, dx1, dx2, dy1, dy2, FoundHits); const Int_t nFound = FoundHits.GetEntries(); if(DebugLevel() > 5) cout << "Row : " << nextRow << " - Number of matching hits found : " << nFound << endl; if(nFound == 0) { if(DebugLevel() > 4) cout << "No hit in row " << nextRow << " for track no. " << i << endl; // check to see if we want to throw track away // When to throw track away : // 1) track has only 1 hit and missed the next row // 2) track cannot get the minimum number of hits const Int_t nRow = fParams->GetNumberOfActiveRows(); if(track->GetNhit() == 1 || (track->GetNhit()+irow-2) < (nRow - fMaxRowsMissed)) { // irow-2 is the maximal number of hits you can get in the // last irow-1 rows. track->MarkAsBad(); if(DebugLevel() > 2) cout << "Track : " << i << " marked as bad" << endl; } } else { // Update track information. // // Split local track when more than one hit is found // and add new hit // The loop over hits has to count down to zero for (Int_t ihit = nFound-1; ihit >= 0; ihit--) { BrTpcTrackCandidate *nextTrack; if (ihit > 0) { // make a copy of old track nextTrack = BrTpcTrackModule::AddNewTrackCandidate(track); // check that track limit is not exceeded if(fStatus) return; if(DebugLevel() > 0 && ihit == 1) { cout << "Track No. " << i << " forked " << nFound << " times - hits on track : " << nextTrack->GetNhit() << " + 1" << endl; } } else { nextTrack = track; } BrTpcHit* hit = (BrTpcHit*)FoundHits.At(ihit); nextTrack->AddHit(hit); // Update position and direction information nextTrack->Fit(); nextTrack->SetLastRow(nextRow); nextTrack->SetStatus(BrTpcTrackCandidate::kOk); } } } // local track loop } // row loop // Remove the bad tracks DeleteBadTracks(); // cout << endl << "Number of track candidates : " << GetNTrackCandidates() // << endl << endl; } //_________________________________________________________________________ void BrTpcTrackFollowModule::FindHitsNextRow(const Int_t rowNo, const Float_t *proj, const Float_t dx1, const Float_t dx2, const Float_t dy1, const Float_t dy2, TObjArray &FoundHits) { // Loop over the hits in row rownumber-1 and accepts the hits in the // interval // proj[0] - dx1 < x < proj[0] + dx2 && // proj[1] - dy1 < y < proj[1] + dy2 // The hits are then stored in the array FoundHits // A lot of diagnosic histograms are filled! if(DebugLevel()>10) { cout << "Cuts x :" << setw(10) << proj[0] << " -" << setw(10) << dx1 << " +" << setw(10) << dx2 <<endl; cout << "Cuts y :" << setw(10) << proj[1] << " -" << setw(10) << dy1 << " +" << setw(10) << dy2 <<endl; } FoundHits.Clear(); const Int_t numHits = fHitsInRow[rowNo-1].GetEntriesFast(); for(Int_t i = 0; i < numHits; i++) { BrTpcHit *hit = (BrTpcHit*)fHitsInRow[rowNo-1].At(i); // hit->Print(); cout << endl; Float_t dx = hit->GetPos()[0] - proj[0]; Float_t dy = hit->GetPos()[1] - proj[1]; if(HistOn()){ hHitSearchDx->Fill(dx); hHitSearchDy->Fill(dy); hHitRowSearchDx->Fill(dx, rowNo); hHitRowSearchDy->Fill(dy, rowNo); } if(-dx1 < dx && dx < dx2 && -dy1 < dy && dy < dy2) { // cout << "Accepted!" << endl; FoundHits.Add(hit); if(HistOn()) { hHitAcceptedDx->Fill(dx); hHitAcceptedDy->Fill(dy); hHitRowAcceptedDx->Fill(dx, rowNo); hHitRowAcceptedDy->Fill(dy, rowNo); } } } } //____________________________________________________________________ void BrTpcTrackFollowModule::Print(Option_t* option) const { // Print module information // See BrModule::Print for options. // In addition this module defines the Option: // <fill in here> BrModule::Print(option); TString opt(option); opt.ToLower(); if (opt.Contains("d")) cout << endl << " Original author: Peter H.L. Christiansen" << endl << " Last Modifications: " << endl << " $Author: cholm $" << endl << " $Date: 2002/01/03 19:53:42 $" << endl << " $Revision: 1.5 $ " << endl << endl << "-------------------------------------------------" << endl; cout << GetName() << "TpcTrackFollowModule" << endl; cout << " SearchWidthX = " << setw(8) << fSearchWidthX << endl; cout << " SearchWidthY = " << setw(8) << fSearchWidthY << endl; cout << " MaxRowsMissed = " << setw(8) << fMaxRowsMissed << endl; BrTpcTrackModule::Print(option); } //____________________________________________________________________ // // $Log: BrTpcTrackFollowModule.cxx,v $ // Revision 1.5 2002/01/03 19:53:42 cholm // Prepared to use BrTableNames class (or perhaps BrDetectorList) for table names // // Revision 1.4 2001/11/02 13:54:05 pchristi // Added changes to BrTpcTrackFollowModule that reflects the use of the // new class BrTrackResidual. // // Revision 1.3 2001/08/14 19:29:33 pchristi // Changed the Print methods to pass the option to BrModule. // Also removed debug statement in BrTpcTrackFollowModule. // // Revision 1.2 2001/07/06 10:24:32 pchristi // Changed assert to the ROOT Assert in TError.h // // Revision 1.1.1.1 2001/06/21 14:55:13 hagel // Initial revision of brat2 // // Revision 1.1 2001/06/17 17:41:44 pchristi // The new tracking classes. // // |
||||||
|
This page automatically generated by script docBrat by Christian Holm |
Copyright ; 2002 BRAHMS Collaboration
<brahmlib@rcf.rhic.bnl.gov>
|