|
//_____________________________________________________________________ // // Drift chamber view combiner // The BrDCViewCombiner class finds local tracks from DC's // //_____________________________________________________________________ // // $Id: BrDCViewCombiner.cxx,v 1.7 2002/02/20 13:53:15 ufstasze Exp $ // $Author: ufstasze $ // $Date: 2002/02/20 13:53:15 $ // // #ifndef ROOT_TString #include "TString.h" #endif #ifndef ROOT_TDirectory #include "TDirectory.h" #endif // #ifndef BRAT_BrDCHitPosition #include "BrDCHitPosition.h" #endif #ifndef BRAT_BrDCHit #include "BrDCHit.h" #endif #include "BrDCViewCombiner.h" #include "BrLocalTrack.h" #include "BrTrackHit2D.h" #include "BrMath.h" ClassImp(BrDCViewCombiner); //_____________________________________________________________ BrDCViewCombiner::BrDCViewCombiner() { SetState(kSetup); // DO NOT USE //fTrackHits = 0; //fTrackHits2D = 0; } //_____________________________________________________________ BrDCViewCombiner::BrDCViewCombiner(const Char_t *Name, const Char_t *Title, BrDCClusterFinder *clusterFinder, const Float_t asigma, const Int_t max_miss_main, const Int_t max_miss_sub, const Int_t max_miss_other, const Float_t uvdx, const Float_t uvdy, Bool_t EventByEvent) : BrModule(Name,Title) { SetState(kSetup); fClusterFinder = clusterFinder; fDCCutMaxMissMain = max_miss_main; fDCCutMaxMissSub = max_miss_sub; fDCCutMaxMissOther = max_miss_other; fUVDx = uvdx; fUVDy = uvdy; fAsigma = asigma; fEventByEvent=EventByEvent; } //_____________________________________________________________ BrDCViewCombiner::~BrDCViewCombiner(){ } //_____________________________________________________________ void BrDCViewCombiner::Init() { SetState(kInit); fSelectedDCHits = new BrDataTable("hits", "hits"); if(!strcmp("T3",GetName()))fXYNumPlane=6; else if (!strcmp("T4",GetName()) || !strcmp("T5",GetName()))fXYNumPlane=4; else cout<<"BrDC2DTrackFinder::Init(): check your BrDCTrackingModule constructor"<<endl; // specific options for matching hits to track. // If fCheckSingleHits is TRUE we look for single // (one plane) hit but only in the case when we coudn't // find the combined hit. // It was found that for T3, this option significantly // slows down the program execution. if(!strcmp("T4",GetName()))fCheckSingleHits = kTRUE; if(!strcmp("T5",GetName()))fCheckSingleHits = kTRUE; if(fBfsAngle>5 || fRunNo>6166){ if(!strcmp("T3",GetName()))fCheckSingleHits = kTRUE; if(!strcmp("T3",GetName()))fStrictParam=1.0; if(!strcmp("T4",GetName()))fStrictParam=0.8; if(!strcmp("T5",GetName()))fStrictParam=0.8; }else{ if(!strcmp("T3",GetName()))fCheckSingleHits = kFALSE; if(!strcmp("T3",GetName()))fStrictParam=2.0; if(!strcmp("T4",GetName()))fStrictParam=0.8; if(!strcmp("T5",GetName()))fStrictParam=0.8; } } //_____________________________________________________________ void BrDCViewCombiner::DefineHistograms() { if (GetState() != kInit) { Stop("DefineHistograms", "Must be called after Init"); return; } // The init member function/method is used to define histograms // The histograms is added to the list of hist objects in this // module. // For combining performance TDirectory* saveDir = gDirectory; TDirectory* histDir = gDirectory->mkdir(Form("%sCombine", GetName())); histDir->cd(); Char_t HistName[64],HistTitle[64]; sprintf(HistName,"CombInU%s",GetName()); sprintf(HistTitle,"XY Track - U Deviation in %s",GetName()); fCombUDuc = new TH1F(HistName,HistTitle,1000,-1.0,1.0); sprintf(HistName,"CombInV%s",GetName()); sprintf(HistTitle,"XY Track - V Deviation in %s",GetName()); fCombVDuc = new TH1F(HistName,HistTitle,1000,-1.0,1.0); gDirectory = saveDir; } //_____________________________________________________________ void BrDCViewCombiner::CombineViews(BrDCViewList *viewList, TObjArray *trackHits2D, TObjArray *trackHits, BrClonesArray *LocalTracks, Int_t &nLocalTracks) { //Combine views // this is only for debugging purpose: //ofstream OUTTr; //Char_t file[64]; //sprintf(file,"./Visu/%stracks.loc",GetName()); //OUTTr.open(file,ios::out); //OUTTr.setf(ios::showpoint); //------------------------------------ Int_t flg; Int_t i,j,k,it,nok; Int_t maxmiss[3]; Int_t NumHitCMBMod; Int_t NumTra2DMod1; BrDataTable Tra2DMod1; Int_t NumTra2DMod2; BrDataTable Tra2DMod2; Float_t u1,u2,u3,du1,du2,du3,z1,z2,z3,au1,au2,th1,th2,th3; Float_t uc,duc,uc1,uc2,x,y,dx2,dy2; Float_t arg1,arg2,arg3; Float_t loctra_pos[3],loctra_vec[3]; Int_t nloctra,numloctra; BrLocalTrack *loctra_p; BrLocalTrack *loctra_p1; Int_t IDSave,IDCurrent,IDCur,IDCurrentPl,nmmiss; Int_t oview; BrDCHit *hit_p; static const Double_t rad = .01745329; // main_view = GetMainView(); // sub_view = GetSubView(); fTrackHits2D = trackHits2D; fTrackHits = trackHits; maxmiss[0] = fDCCutMaxMissMain; maxmiss[1] = fDCCutMaxMissSub; maxmiss[2] = fDCCutMaxMissOther; flg = 0; GetSelected2DTracks(0,flg,&Tra2DMod1); flg = 1; GetSelected2DTracks(0,flg,&Tra2DMod2); //if(DebugLevel()>5) cout<<"Getting number of trahitmods"<<endl; NumTra2DMod1 = Tra2DMod1.GetEntries(); NumTra2DMod2 = Tra2DMod2.GetEntries(); if(fEventByEvent){ cout<<"Number of 2D X Tracks: "<<NumTra2DMod1<<endl; cout<<"Number of 2D Y Tracks: "<<NumTra2DMod2<<endl; } //Int_t ii;cin>>ii; for(it=0;it<NumTra2DMod1;it++) { BrTrackHit2D* tra2dd = (BrTrackHit2D*)Tra2DMod1.At(it); u1 = tra2dd->GetPos()[1]; du1 = tra2dd->GetDpos(); z1 = tra2dd->GetPos()[0]; au1 = tra2dd->GetVec(); th1 = (Float_t)viewList->GetMainView(); for(j=0;j<NumTra2DMod2;j++) { BrTrackHit2D* tra2dd2 = (BrTrackHit2D*)Tra2DMod2.At(j); u2 = tra2dd2->GetPos()[1]; du2 = tra2dd2->GetDpos(); z2 = tra2dd2->GetPos()[0]; au2 = tra2dd2->GetVec(); th2 = (Float_t)viewList->GetSubView(); nloctra = LocalTracks->GetLast() + 1; BrClonesArray &loctra = *LocalTracks; loctra_p = new(loctra[nLocalTracks++]) BrLocalTrack(); loctra_pos[0] = u1; loctra_pos[1] = u2 + au2*(z1-z2); loctra_pos[2] = z1; loctra_p->SetPos(loctra_pos); loctra_vec[0] = au1; loctra_vec[1] = au2; loctra_vec[2] = (Float_t)1.0; loctra_p->SetVec(loctra_vec); loctra_p->SetNhit(tra2dd->GetNhit() + tra2dd2->GetNhit()); loctra_p->SetFlg(0); loctra_p->SetTr1(tra2dd); loctra_p->SetTr2(tra2dd2); loctra_p->SetID((Int_t)loctra_p); //cout<<au1<<" "<<au2<<" "<<u1-au1*z1<<" "<<u2-au2*z2<<endl; // OUTTr<<loctra_p->GetVec()[0]<<" "<<loctra_p->GetVec()[1]<<" " // <<loctra_p->GetPos()[0]-au1*z1<<" "<<loctra_p->GetPos()[1]-au2*z1<<endl; numloctra = nloctra; IDSave = numloctra; IDCurrent = numloctra; nmmiss = 0; //Loop over intermediate planes u, v for verification Int_t nov = viewList->GetNumPlaneViews(2); //2 to indicate other view //cout<<it<<" "<<j<<" viewList->GetNumPlaneViews(2): "<<nov<<endl; for(oview=0;oview<nov;oview++) { IDCurrentPl = IDCurrent; nok = 0; th3 = viewList->GetViewPlaneAngle(2,oview+1); //2 is to indicate other view TObjArray *detectorHitsOtherView = viewList->GetDetectorHits(2,oview+1); NumHitCMBMod = detectorHitsOtherView->GetEntries(); //cout<<"NumHitCMBMod= "<<NumHitCMBMod<<endl; for(k=0;k<NumHitCMBMod;k++) { hit_p = (BrDCHit*)detectorHitsOtherView->At(k); u3 = hit_p->GetPos()[0]; z3 = hit_p->GetPos()[2]; du3 = hit_p->GetDpos()[0]; // cout<<"In CombineViews: "<<oview<<" "<<th1<<" "<<th2<<" "<<th3<<endl; uc1 = u1 + au1*(z3-z1); //normally x-view uc2 = u2 + au2*(z3-z2); //normally y-view uc = uc1*(Float_t)TMath::Cos(th3*rad) - uc2*(Float_t)TMath::Sin(th3*rad); //duc=fUVDx; //duc=du3 + 0.5*(du2+du1)*0.309/(Float_t)1.41; //.309=sin(18) (default) duc=du3 + (du2+du1)*0.309/(Float_t)1.41; //.309=sin(18) (T3 test) //duc=du3; // to have previous version //if(!strcmp("T5",GetName())) //cout<<"Comb. View: "<<arg1<<" "<<arg2<<" "<<arg3<<" "<<duc<<endl; /* if(DebugLevel()>4) { cout<<"uc,duc,uc-sigma*duc,u3,uc+sigma*duc "; cout<<uc<<","<<duc<<","; cout<<uc-fAsigma*duc<<","<<u3<<","<<uc+fAsigma*duc<<endl; } */ if(th3==-18)fCombUDuc->Fill(uc - u3); if(th3==18) fCombVDuc->Fill(uc - u3); Float_t lowerLimit = uc-fAsigma*duc; Float_t upperLimit = uc+fAsigma*duc; //cout<<"nok,NMiss "<<nok<<" "<<nmmiss<<" "<<lowerLimit<<","<<u3<<","<<upperLimit<<endl; //if(BrMath::BTWN(uc-fAsigma*duc,u3,uc+fAsigma*duc)) { if(BrMath::BTWN(lowerLimit,u3,upperLimit)) { nok++; if(nok == 1) { // Not the right place to update the #hits on this // track?? // // The loctra_ID need not be the same throughout the // whole loop. // Update the #hits information and save the current // hits for all the local tracks for(i=IDSave;i<=IDCurrent;i++) { BrLocalTrack *loctra_p2 = (BrLocalTrack*)LocalTracks->At(i); loctra_p2->IncNhit(); InsertTrahit(hit_p,loctra_p2); } } //else if(nok > 1) { else if(nok == 1) { //if(fEventByEvent)cout<<"New Ghost Tracks: "<<IDCurrentPl-IDSave+1<<endl; //Add ghost tracks for(i=IDSave;i<=IDCurrentPl;i++) { //loctra_p1 = new BrLocalTrack(loctra_p); //LocalTracks->Add(loctra_p1); BrClonesArray &loctra = *LocalTracks; // original version loctra_p1 = new(loctra[nLocalTracks++]) BrLocalTrack(loctra_p); loctra_p1->SetID((Int_t)loctra_p1); BrLocalTrack* loctra_pi = (BrLocalTrack*)LocalTracks->At(i); CopyTrahit1(loctra_pi,loctra_p1);//??? } IDCur = IDCurrent; IDCurrent = LocalTracks->GetLast(); for(i=IDCur+1;i<=IDCurrent;i++) { BrLocalTrack *loctra_p2 = (BrLocalTrack*)LocalTracks->At(i); InsertTrahit(hit_p,loctra_p2); } } } } //loop over k //////////////////////////////////////////////////////////// for(Int_t ipv=0;ipv<2;ipv++){ if(nok==0 && fCheckSingleHits){ //try to find something in single planes Int_t idet = oview/2; Int_t iview = oview%2 + 2; Int_t ipl = (fXYNumPlane-1) + (iview-1)*2 + ipv; z3 = (fClusterFinder->GetDCDriver())->GetPlanePosition(idet,ipl); uc1 = u1 + au1*(z3-z1); //normally x-view uc2 = u2 + au2*(z3-z2); //normally y-view uc = uc1*(Float_t)TMath::Cos(th3*rad) - uc2*(Float_t)TMath::Sin(th3*rad); fClusterFinder->SelectDCHits(idet,ipl,fSelectedDCHits); Int_t numSelectedHits = fSelectedDCHits->GetEntries(); for(Int_t k=0;k<numSelectedHits;k++) { BrDCHitPosition *hitpos_p = (BrDCHitPosition*)fSelectedDCHits->At(k); BrDcDigInfo *diginfo_p = hitpos_p->GetDigDCInfo(); if(diginfo_p->GetUsed() == 0) { u3 = hitpos_p->GetUhit(); du3 = (fClusterFinder->GetDCParams())->GetPosres(); duc = du3/fStrictParam; // to be rather strict (was 2.9) Float_t lowerLimit = uc-fAsigma*duc; Float_t upperLimit = uc+fAsigma*duc; if(BrMath::BTWN(lowerLimit,u3,upperLimit)) { //Coming here means that we found a hit that is consistent //with this 2D track, so first create ne BrDCHit: BrClonesArray &hit = *fClusterFinder->GetDetectorHits(); Int_t nhit = (fClusterFinder->GetDetectorHits())->GetLast()+1; hit_p = new(hit[nhit]) BrDCHit(); Float_t hit_pos[3]; Float_t hit_dpos[3]; hit_pos[0] = u3; hit_pos[1] = (Float_t)0.0; hit_pos[2] = z3; hit_p->SetPos(hit_pos); hit_dpos[0] = du3; hit_dpos[1] = (Float_t)0.0; hit_p->SetDpos(hit_dpos); //cout<<"VC: "<<idet<<" "<<ipl<<" "<<iview<<" "<<oview<<endl; hit_p->SetImod((idet + 1)*100 + iview+1); hit_p->SetStat(1); // save primary hit hit_p -> SetNumCombined(1); hit_p -> SetU(u3,0,0); hit_p -> SetZ(z3,0,0); hit_p -> SetPlane(ipl,0,0); hit_p -> SetWire(hitpos_p->GetNwire(),0,0); hit_p -> SetTdc(hitpos_p->GetTdc(),0,0); nok++; if(nok == 1) { for(i=IDSave;i<=IDCurrent;i++) { BrLocalTrack *loctra_p2 = (BrLocalTrack*)LocalTracks->At(i); loctra_p2->IncNhit(); InsertTrahit(hit_p,loctra_p2); } diginfo_p->IncUsed(); } } //if(BrMath::BTWN(lowerLimit,u3,upperLimit)) } } // loop over k fSelectedDCHits->Clear(); } // if(nok==0 && fSingleHitsCheckingOption) } // loop over ipv /////////////////////////////////////////////////////////// if(nok == 0) { nmmiss++; if(nmmiss>maxmiss[2]) { for(i=IDSave;i<=IDCurrent;i++) { loctra_p = (BrLocalTrack*)LocalTracks->At(i); //RemoveTrack2(loctra_p); loctra_p->MarkAsBad(); } break; } } //detectorHitsOtherView->Clear(); } //loop over oview // if(DebugLevel() > 2) cout<<"nmiss for loctra = "<<nmmiss<<endl; } //loop over j (y view tracks) } //loop over i (x view tracks) //set track status to be 1 //if(DebugLevel()>1) cout<<"NumLoctra before cleaning is "<<LocalTracks->GetEntries()<<endl; DeleteBadTrack(LocalTracks,nLocalTracks); // if(DebugLevel()>1) cout<<"NumLoctra after cleaning is "<<LocalTracks->GetEntries()<<endl; //Clear container arrays Tra2DMod1.Clear(); Tra2DMod2.Clear(); //Tra2DMod1.Delete(); //Tra2DMod2.Delete(); //OUTTr.close(); } void BrDCViewCombiner::GetSelected2DTracks(const Int_t stat, const Int_t flg, BrDataTable *Tra2DMod) { //Get selected 2DTracks and put into BrDataTable //stat is status of 2D track to be selected //flg is view (main or sub; usually x or y) Int_t i; Int_t NumTra2D; BrTrackHit2D *tra2d_p; Tra2DMod->Clear(); //Tra2DMod->Delete(); NumTra2D = fTrackHits2D->GetLast()+1; for(i=0;i<NumTra2D;i++) { tra2d_p = (BrTrackHit2D*)fTrackHits2D->At(i); if(( tra2d_p->GetStat() == stat || stat == IANY ) && ( tra2d_p->GetFlg() == flg || flg == IANY )) { Tra2DMod->Add(tra2d_p); } } } //_____________________________________________________________________________ void BrDCViewCombiner::InsertTrahit(BrDCHit *hit, BrLocalTrack* localtrack) { // Fill a new entry in the mapping 'table' TrackHits with the // pointers to a detector hit (hit) and to a localtrack. // BrTrackHit* trahit_p; trahit_p = new BrTrackHit(); fTrackHits->Add(trahit_p); trahit_p->SetImod(hit->GetImod()); trahit_p->SetTrackID(localtrack->GetID()); trahit_p->SetHitID(hit->GetID()); trahit_p->SetLoctra(localtrack); trahit_p->SetDetectorHit(hit); localtrack->AddTrackHit(trahit_p); } //_____________________________________________________________ void BrDCViewCombiner::CopyTrahit1(BrLocalTrack* OldTr,BrLocalTrack* NewTr) { //Copy trahit entries excluding the last entry in hitcmb. Int_t itr; Int_t NumTraHitFound; TObjArray TraHitFound; BrTrackHit* trahit_p; BrDCHit* hit_p; Int_t NumTraHit; NumTraHitFound = 0; //First, figure out where the last entry is NumTraHit = fTrackHits->GetLast()+1; for(itr=0;itr<NumTraHit;itr++) { trahit_p = (BrTrackHit*)fTrackHits->At(itr); if( trahit_p->GetLocalTrack() == OldTr ) { TraHitFound.Add(trahit_p); } } NumTraHitFound = TraHitFound.GetEntries(); for(itr=0;itr<NumTraHitFound-1;itr++) { //hit_p = ((BrTrackHit*)TraHitFound.At(itr))->GetDetectorHit(); hit_p = (BrDCHit*)((BrTrackHit*)TraHitFound.At(itr))->GetDetectorHit(); if(hit_p != 0) { /* BrClonesArray &trahit = *TrackHits; Int_t ntrahit = TrackHits->GetLast()+1; trahit_p = new(trahit[ntrahit]) BrTrackHit(); */ trahit_p = new BrTrackHit(); fTrackHits->Add(trahit_p); //trahit_p->SetHitcmb(hitcmb_p); trahit_p->SetDetectorHit(hit_p); trahit_p->SetTrackID(NewTr->GetID()); trahit_p->SetLoctra(NewTr); // Set up relation with loctra NewTr->AddTrackHit(trahit_p); } } TraHitFound.Clear(); //TraHitFound.Delete(); } //_____________________________________________________________________________ void BrDCViewCombiner::RemoveTrack(BrLocalTrack* localtrack) { //For future reference, this can be done in a much easier way //using the numtrahit and trahit entries in loctra. That gives //us the information directly without having to do a search. Int_t i; BrTrackHit *trahit_p; BrDCHit *hit_p; Int_t NumTraHit; NumTraHit = fTrackHits->GetEntries(); for( i=0;i<NumTraHit;i++) { trahit_p = (BrTrackHit*)fTrackHits->At(i); if(trahit_p->GetLocalTrack() == localtrack) { // First, find if any hitcmb corresponds to this and reset stat hit_p = (BrDCHit*)trahit_p->GetDetectorHit(); hit_p->SetStat((hit_p->GetStat())%16 + (hit_p->GetStat()/16-1)*16); localtrack->RemoveTrackHit(trahit_p); fTrackHits->RemoveAt(i); } } fTrackHits->Compress(); localtrack->MarkAsBad(); } //_____________________________________________________________________________ void BrDCViewCombiner::DeleteBadTrack(BrClonesArray *localTracks, Int_t &nLocalTracks) { Int_t i; Int_t NumLoctra; BrLocalTrack *loctra_p; NumLoctra = localTracks->GetEntries(); for(i=0;i<NumLoctra;i++) { loctra_p = (BrLocalTrack*)localTracks->At(i); if( loctra_p->IsBad()) { localTracks->RemoveAt(i); nLocalTracks--; } } localTracks->Compress(); } // $Log: BrDCViewCombiner.cxx,v $ // Revision 1.7 2002/02/20 13:53:15 ufstasze // fCheckSingleHits and fStrictParam set according to spectr. angle and run number. // // Revision 1.6 2002/02/15 16:05:56 ufstasze // fCheckSingleTrack set to kTRUE // // Revision 1.5 2001/11/12 10:24:52 ufstasze // new values for fStrictParam for T4 and T5 // // Revision 1.4 2001/08/23 15:53:45 staszel // Now, we don't use RemoveTracks method. We just // mark the track as bad. Is significantly saves cpu and // decreases memory leak. // // Revision 1.3 2001/08/23 10:52:32 staszel // CheckSingleHits option has been introduced. // // Revision 1.2 2001/08/07 13:48:49 ouerdane // Cleaned up a lot the DC tracking modules and utilities: // Made all utilities deriving from BrModule // Set the class version number to 0 (NOT 1 like before) // Removed + in the LinkDef file for all modules // Added ClassDef in the View Combiner class (didn't exist) // Structurized the histogramming section in TDirectories // // Revision 1.1.1.1 2001/06/21 14:55:12 hagel // Initial revision of brat2 // // Revision 1.7 2001/06/20 12:11:45 staszel // BrDetectorHit replaced by BrDCHit. // // Revision 1.6 2001/06/13 13:20:39 staszel // father development of dc tracking // // Revision 1.5 2001/01/19 17:55:06 staszel // Changes having to do with development of BrDCTrackingModule // // Revision 1.4 2000/10/25 17:58:06 hagel // Further development of BrDCTrackingModule // // Revision 1.3 2000/10/05 02:21:31 hagel // Further development of BrDCTrackingModule // // Revision 1.2 2000/09/29 02:15:59 hagel // Changes having to do with development of BrDCTrackingModule // // Revision 1.1 2000/09/19 21:10:55 hagel // Initial placeholder revision // |
||||||
|
This page automatically generated by script docBrat by Christian Holm |
Copyright ; 2002 BRAHMS Collaboration
<brahmlib@rcf.rhic.bnl.gov>
|