|
//____________________________________________________________________________ // // BrMrsTrackingModule // // Tracking module for the MRS. The module will setup the module that combines // the local tracks found front and after the D5 magnet. // // Requires: BrDetectorTracks Table in inputnode of the // two tracking TPS TPM1 and TPM2. // // Note that this new version always get the fields from the database // Please use it being aware of this. // // This may be a real issue when running simulated data - but someone will discover this // one day - be forwarned. (fv) //____________________________________________________________________________ // // $Id: BrMrsTrackingModule.cxx,v 1.10 2002/08/30 17:36:34 videbaek Exp $ // $Author: videbaek $ // $Date: 2002/08/30 17:36:34 $ // #include <BrIostream.h> #include <iomanip> #include <BrMrsTrackingModule.h> #include <TMath.h> #include <TDirectory.h> #include <TString.h> #include "BrVector3D.h" #include "BrLine3D.h" #include "BrPlane3D.h" #include "BrEvent.h" #include "BrModuleMatchTrack.h" #include "BrDetectorTrack.h" #include "BrGeometryDbManager.h" #include "BrParameterDbManager.h" #include "BrRunInfoManager.h" //______________________________________________________________________ ClassImp(BrMrsTrackingModule); //______________________________________________________________________ BrMrsTrackingModule::BrMrsTrackingModule() { //// default constructor SetState(kSetup); fCombineModule = 0; SetUseDb(); } //______________________________________________________________________ BrMrsTrackingModule::BrMrsTrackingModule(const Char_t* Name, const Char_t* Title) : BrModule(Name, Title) { // // Normal constructor. // This Name & title is informational. // SetState(kSetup); fCombineModule = new BrModuleMatchTrack("MRS","MRS matching", "TPM1","TPM2","D5"); fCombineModule->SetHistMomentumRange(3.0); SetUseDb(); } //______________________________________________________________________ BrMrsTrackingModule::~BrMrsTrackingModule() { //// default destructor if(fCombineModule) delete fCombineModule; } //__________________________________________________________________________ void BrMrsTrackingModule::Init(){ // // Clear datastructures // request geometry information using the GeometryDbManager. // SetState(kInit); fCombineModule->SetDebugLevel(DebugLevel()); fCombineModule->SetVerbose(Verbose()); fCombineModule->Init(); if(fCombineModule->GetStatus()== BrModule::kAbort){ Abort("Init","Combine Module Failed"); return; } fFrontVolume = fCombineModule->GetFrontVolume(); fBackVolume = fCombineModule->GetBackVolume(); fMagnetVolume = fCombineModule->GetMagnetVolume(); BrGeometryDbManager *geomDb = BrGeometryDbManager::Instance(); BrParameterDbManager *parDb = BrParameterDbManager::Instance(); fTofParams = (BrDetectorParamsTof*)parDb-> GetDetectorParameters("BrDetectorParamsTof", "TOFW"); fPanelPar = new BrDetectorParamsTof* [fTofParams->GetNoPanels()]; fPanelVol = new BrDetectorVolume* [fTofParams->GetNoPanels()]; for (Int_t p = 0; p < fTofParams->GetNoPanels(); p++) { fPanelPar[p] = (BrDetectorParamsTof*)parDb-> GetDetectorParameters("BrDetectorParamsTof", Form("TFP%d", p+1)); fPanelVol[p] = (BrDetectorVolume*)geomDb-> GetDetectorVolume("BrDetectorVolume", Form("TFP%d", p+1)); } } //__________________________________________________________________________ void BrMrsTrackingModule::DefineHistograms() { // // define histograms for modules owned by this 'super module' // TDirectory* savdir = gDirectory; TDirectory* hdir = savdir->mkdir("Tracking_MRS"); if (!hdir) { Warning("DefineHistograms", "Could not create histogram subdirectory"); return; } hdir->cd(); fCombineModule->Book(); hGhost = new TH2F("hGhost", "Tracks after ghostbusting vs Tracks before", 6, 0.5, 6.5, 6, 0.5, 6.5); hStatus = new TH1F("hStatus", "Ghost busting status", 10, 0.5, 10.5); hTargetTz = new TH1F("hTz","MRS target image (z)", 200, -50., 50.); hTargetTx = new TH1F("hTx","MRS target image (x)", 200, -10., 10.); hTargetTy = new TH1F("hTy","MRS target image (y)", 200, -10., 10.); hP = new TH1F("hMomentum", "Accepted MRS tracks momentum",320, -4., 4.); hTheta = new TH1F("hTheta","MRS Track Theta", 200, 0.0, 2.0); hPhi = new TH1F("hPhi","MRS Track Phi", 200, -.2, .2); hPathLength = new TH1F("hPathLength","MRS path length",400,300.,500.); hPartialPath = new TH1F("hPartialPath","MRS path length from TPM1 to TOFW", 400, 100., 500.); hSlatPos = new TH1F("hSlatPos","Slat no for projected tracks", 130,-0.5, 129.5); hTrkBbDist = new TH1F("hTrkDistBb", "Distance between track projection to BB " "plane normal to MRS", 400, -10, 10); hTrkBbXY = new TH2F("hTrkBbXY", "Track projection to BB vertex (small tubes)" "plane normal to MRS (in local coord.)", 200, -10, 10, 200, -10, 10); hTrkBbZ = new TH1F("hTrkBbZ", "BB Vertex (small tubes) - Track intersection with beam line", 500, -25, 25); // Restore directory savdir->cd(); } //___________________________________________________________________ void BrMrsTrackingModule::SetMagnetField() { // ---------------------------------------------------- // get field information from database // update info each time a field is needed (event loop) // ---------------------------------------------------- BrRunInfoManager* runMan = BrRunInfoManager::Instance(); if (!runMan) { Abort("SetMagnetField", "Could not get an instance of BrRunInfoManager"); return; } const BrRunInfo* run = runMan->GetCurrentRun(); if(!run){ Abort("SetMagnetField", "Could not get current run"); return; } if (run->GetRunNo() == -1) { Abort("SetMagnetField", "Run number is -1, check it out"); return; } // get magnet setting fMagnetVolume->SetCurrent(run->GetD5Set(), run->GetD5Pol()); if (Verbose() > 10) cout << " *** BrMrsTrackingModule::SetMagnetField() -- > D5: " << fMagnetVolume->GetField() << endl; } //___________________________________________________________________________ void BrMrsTrackingModule::Begin() { // Set magnet field and pick up the platform angle // this is done at begin time since one might scan sequences // from diff. runs with diff. settings (maybe not a good thing to do) if(fUseDb) SetMagnetField(); } //___________________________________________________________________________ void BrMrsTrackingModule::Event(BrEventNode* InputNode, BrEventNode* OutputNode) { // // Normal Event member function for MRS tracking // // Requires: Two instances of BrDataTable named "DetectorTrack TPM1" and // "DetectorTrack TPM2" in either InputNode or OutputNode. // These tables have to be filled with instances of BrDetectorTrack // from their respective TPCs. // Result: BrDataTable "MrsTracks" added to outputnode. // The elements of this table are BrMrsTrack. // // The module has the following steps: // a) Combine Front (TPM1) and Back (TPM2) tracks // b) Track combined tracks back to nominal vertex // c) Check goodness of tracks and perform cleanup // of ghost tracks, i.e. cases where the same track // are used several times. The tracks are not removed // but are marked as 'ghosts'. // d) Add tracks to OutputNode (good tracks and potential // ghosts if requested). // SetState(kEvent); SetStatus(kOk); fCombineModule->Clear(); // -------------------------------------------- // -- Combining local tracks -- // -------------------------------------------- BrEventNode* node = InputNode; // YACK!! if (!InputNode->GetDataTable("DetectorTrack TPM1") || !InputNode->GetDataTable("DetectorTrack TPM2")) node = OutputNode; fCombineModule->CombineFrontBack(node); if(DebugLevel() > 2) cout << "<BrMrsTrackingModule::Event>: Get array of matched tracks" << endl; // ---------------------------------------- // -- Get matched tracks -- // ---------------------------------------- fTracksFrontBack = fCombineModule->GetMatchedTracks(); if(DebugLevel() > 1) cout << "<BrMrsTrackingModule::Event>: Number of matched tracks = " << fTracksFrontBack->GetEntries() << endl; if (!fTracksFrontBack->GetEntries()) return; // ---------------------------------------- // -- Prepare output table -- // ---------------------------------------- BrDataTable* mrstable = new BrDataTable("MrsTracks"); OutputNode->AddObject(mrstable); // ---------------------------------------- // -- Track down to primary vertex -- // ---------------------------------------- TrackToVertex(mrstable); // ----------------------------------------------------------------------------- // -- Pick up BB vertex (the new one) if required -- // -- (Note the increasing usage of the check both in input and output node -- // -- for data objects. This shows the wrong usage of the event model -- // -- and bratmain. No good.) -- // ----------------------------------------------------------------------------- if (HistOn()) { BrBbVertex* bbVtx = (BrBbVertex*)InputNode->GetObject("BB Vertex"); if (!bbVtx) bbVtx = (BrBbVertex*)OutputNode->GetObject("BB Vertex"); // check only with small tubes vtx // well maybe not for pp !! if (bbVtx) if (bbVtx->GetMethod() == 2) TrackToBbVertex(mrstable, bbVtx); } // ---------------------------------------- // -- Do some selected print out -- // ---------------------------------------- for (Int_t i = 0; i < mrstable->GetEntries(); i++) { if (Verbose() > 25) ((BrMrsTrack*)mrstable->At(i))->Print("A"); else if (Verbose() > 20) ((BrMrsTrack*)mrstable->At(i))->Print("PD"); else if (Verbose() > 15) ((BrMrsTrack*)mrstable->At(i))->Print("P"); else if (Verbose() > 10) ((BrMrsTrack*)mrstable->At(i))->Print(""); } } //___________________________________________________________________ void BrMrsTrackingModule::Finish() { // Finisn entry method // Call combine module finish fCombineModule->Finish(); } //___________________________________________________________________ void BrMrsTrackingModule::TrackToVertex(BrDataTable* tab) { // This routine will loop over all tracks found using their position // in T1 as a starting point for backward tracking. // The target Tx ,Ty and Tz are in the IP coordinate system ie at // x = 0 // Calculate scattering angles from vertex and position // in tpm1. This should be considered a first approximation for // vertex , and particle theta, phi. It should be redone when // checked against vertex positions determined by the vertex modules. // It also calculate the partial tracklength from the vertex to TPM1 // position. // The vertex assumed is that of each individual track not any global // vertexing; this should be done later. // Some diagnostics plots are also filled. // The front back tracks are added to table at this point. if(DebugLevel() > 0) cout << " TrackToVertex: " << endl << " #Tracks = " << fTracksFrontBack->GetEntries() << endl; // // Setup the beamline plane use to project each track // to find the tvertex for this track. // BrPlane3D beamPlane(0,0,0, 0,0,1, 0,1,0); const BrCoordinateSystem* tpm1system = fFrontVolume->GetCoordinateSystem(); Int_t nmatch = fTracksFrontBack->GetEntries(); Int_t nGood = 0; for(Int_t it = 0; it < nmatch; it++) { // pick up matched track BrMatchedTrack *matchTrk = (BrMatchedTrack*) fTracksFrontBack->At(it); if (HistOn()) hStatus->Fill(matchTrk->GetStatus()); // check ghostbusting status // BUT do not reject tracks !!! // if (matchTrk->GetStatus() == 1) nGood++; BrDetectorTrack *frTrk = matchTrk->GetFrontTrack(); BrVector3D tpm1position = tpm1system->TransformToMaster(frTrk->GetPos()); BrLine3D frontGlobalTrack(tpm1position, tpm1system->RotateToMaster(frTrk->GetAlpha())); BrVector3D vertex = beamPlane.GetIntersectionWithLine(frontGlobalTrack); BrVector3D trackdir = tpm1position - vertex; Double_t theta = trackdir.Theta(); Double_t phi = trackdir.Phi(); if(HistOn()) { hTargetTx->Fill(vertex(0)); hTargetTy->Fill(vertex(1)); hTargetTz->Fill(vertex(2)); hPhi->Fill(phi); hTheta->Fill(theta); } // // Create a new MRS track BrMrsTrack* mrstrack = new BrMrsTrack(); mrstrack->SetTrackId(it); mrstrack->SetMatchedTrack(matchTrk); mrstrack->SetTrackVertex(vertex); mrstrack->SetTheta(theta); mrstrack->SetPhi(phi); // check now the other end (TOFW stuff) // Similar to previous trackToTofw Float_t length = PartialPath(mrstrack); mrstrack->SetPartialPath(length); mrstrack->SetPathLength(length + trackdir.Norm()); tab->Add(mrstrack); if (HistOn()) { hPathLength->Fill(mrstrack->GetPathLength()); hPartialPath->Fill(mrstrack->GetPartialPath()); hP->Fill(mrstrack->GetMomentum()); } } if (HistOn()) hGhost->Fill(nmatch, nGood); } //___________________________________________________________________ void BrMrsTrackingModule::TrackToBbVertex(BrDataTable* tab, BrBbVertex* bbVtx) { // protected // create plane centered on Z0 and normal to TPM1 front plane // in global coord.!! // (fv) Why! The relevant resolution of the BbVertex is along the // Z axis so differences should be measured along this axis. On the // normal plane the resolution goes as sigma(z)/sin(mrsangle) // BrVector3D normal(fFrontVolume->GetFrontPlane().GetA(), fFrontVolume->GetFrontPlane().GetB(), fFrontVolume->GetFrontPlane().GetC()); BrVector3D beamDir(0, 0, 1); BrVector3D xDir(1, 0, 0); BrVector3D bb(0, 0, bbVtx->GetZ0()); BrVector3D ip(0, 0, 0); // BB vertex in TPM1 front plane coord. BrVector3D bbLocal; fFrontVolume->GlobalToLocal(bb, bbLocal, 0); BrPlane3D bbPlane(bb, normal); // parallel to TPM1 front plane BrPlane3D longPlane(ip, xDir); // sits at X (global) = 0 for (Int_t t = 0; t < tab->GetEntries(); t++) { BrMrsTrack* mrstrack = (BrMrsTrack*)tab->At(t); BrDetectorTrack* trk = mrstrack->GetFrontTrack(); // project track line on BB plane normal to MRS axis BrLine3D trkLine(fFrontVolume->LocalToGlobal(trk->GetTrackLine())); BrVector3D inters, gInters; gInters = bbPlane.GetIntersectionWithLine(trkLine); // calculate distance from intersection to this point: Float_t dist = (bb - gInters).Norm(); hTrkBbDist->Fill(dist); // convert into plane coord. system fFrontVolume-> GlobalToLocal(gInters, inters, 0); hTrkBbXY->Fill(inters(0) - bbLocal(0), inters(1)); // projection on plane normal to global X axis BrVector3D proj = longPlane.GetIntersectionWithLine(trkLine); hTrkBbZ->Fill(bb(2) - proj(2)); } } //___________________________________________________________________ Float_t BrMrsTrackingModule::PartialPath(BrMrsTrack* track) { // Projects to tofw and // calculate the length from the track position in TPM1 // to the time of flight detector. // If the track does not point within any panel A slatNo of 0 isassigned // and the position is undefined. if(DebugLevel() > 0) cout << " In PartialLength() " << endl; // initial stuff BrDetectorTrack* dbtrack = track->GetBackTrack(); BrDetectorTrack* dftrack = track->GetFrontTrack(); Float_t length = 0; Int_t slat = 0; Int_t panel = -1; BrVector3D local; BrVector3D trkProj; // ------------------- for (Int_t panelNo = 0; panelNo < fTofParams->GetNoPanels(); panelNo++) { BrPlane3D frontPanel(fPanelVol[panelNo]->GetFrontPlane()); BrLine3D trkLine(fBackVolume->LocalToGlobal(dbtrack->GetTrackLine())); trkProj = frontPanel.GetIntersectionWithLine(trkLine); fPanelVol[panelNo]->GlobalToLocal(trkProj, local, 0); if(TMath::Abs(local(0)) < fPanelVol[panelNo]->GetSize()[0]/2. && TMath::Abs(local(1)) < fPanelVol[panelNo]->GetSize()[1]/2.) { panel = panelNo; break; } } if (panel < 0) { BrVector3D wrong(-999, -999, -999); track->SetPointedSlat(0); track->SetPointedPanel(panel); track->SetProjOnTof(wrong); return 0; } slat = fPanelPar[panel]->GetSlatNo(local(0)); if (panel > 0) for(Int_t j = 0; j < panel;j ++) slat += fPanelPar[j]->GetNoSlats(); // get matched track length (from TPM1 to TPM2 track pos) Double_t mtLength = track->GetMatchedTrack()->GetMatchedTrackLength(); // back track position in TPM2 (global c.s.) BrVector3D bkPos; fBackVolume->LocalToGlobal(dbtrack->GetPos(), bkPos, 0); // TPM2 to TOFW BrVector3D backTrackToTof(bkPos - trkProj); length = mtLength + backTrackToTof.Norm(); track->SetPointedPanel(panel); track->SetPointedSlat(slat); track->SetProjOnTof(trkProj); // global coord. sys. if (HistOn()) hSlatPos->Fill(slat); return length; } //_________________________________________________________________________ // // $Log: BrMrsTrackingModule.cxx,v $ // Revision 1.10 2002/08/30 17:36:34 videbaek // Corrected the orderof verbose selection 10,15,20,25.. // // Revision 1.9 2002/06/13 00:42:22 videbaek // Change pathlength in histogram // // Revision 1.8 2002/02/26 08:58:20 cholm // Added const'ness to many methods - needed by ISO C++ // // Revision 1.7 2002/01/27 21:15:31 videbaek // Remove references to hTriggers. This was only used by me for // so checking (could and should have been done outside this module). // // Revision 1.6 2002/01/25 17:03:33 videbaek // Take care of errors from combine module // // Revision 1.5 2001/12/07 21:09:22 videbaek // Add SetUseDb() method to control if magnet filed is picked from Db or Not. // default is to do so, but there is a needed particular for MC studies not to use // the rundatabase. This is a lesson for other modules too. // // Added finish method() // // Revision 1.4 2001/12/06 19:53:58 jens // Updated BrMrsTrackingModule::Event to look for the TPC track tables in the output node if not found in the input node. // Also cleaned up this function a little bit. // // Revision 1.3 2001/11/19 15:07:40 ouerdane // removed useless member fMrsAngle // // Revision 1.2 2001/11/08 21:39:47 videbaek // Do not skip adding tracks even if ghost- should not be taken out- // Just count good. Track removal should be in later stage of analysis. // // Revision 1.1 2001/11/05 06:42:21 ouerdane // added new classes replacing older spectrometer track modules // |
||||||
This page automatically generated by script docBrat by Christian Holm |
Copyright ; 2002 BRAHMS Collaboration
<brahmlib@rcf.rhic.bnl.gov>
|