FEDRA emulsion software from the OPERA Collaboration
EdbShowAlg_GS Class Reference

#include <EdbShowAlg_GS.h>

Inheritance diagram for EdbShowAlg_GS:
Collaboration diagram for EdbShowAlg_GS:

Public Member Functions

Double_t CalcIP (EdbSegP *s, double x, double y, double z)
 
Double_t CalcIP (EdbSegP *s, EdbVertex *v)
 
TObjArray * CheckCleanPairs (EdbSegP *InBT, TObjArray *RecoShowerArrayFromFindPairs)
 
Bool_t CheckPairDuplications (Int_t SegPID, Int_t SegID, Int_t Seg2PID, Int_t Seg2ID, TArrayI *SegmentPIDArray, TArrayI *SegmentIDArray, TArrayI *Segment2PIDArray, TArrayI *Segment2IDArray, int RecoShowerArrayN)
 
 ClassDef (EdbShowAlg_GS, 1)
 
void Convert_InVtxArray_To_InBTArray ()
 
void Convert_TracksArray_To_InBTArray ()
 
 EdbShowAlg_GS ()
 
 EdbShowAlg_GS (EdbPVRec *ali)
 
void Execute ()
 
void Finalize ()
 
TObjArray * FindPairs (EdbSegP *InBT, EdbPVRec *eAli_Sub)
 
TObjArray * GetInVtxArray () const
 
Int_t GetInVtxArrayN () const
 
void Init ()
 
void Initialize ()
 
Bool_t IsPossibleFakeDoublet (EdbSegP *s1, EdbSegP *s2)
 
void SetInVtx (EdbVertex *vtx)
 
void SetInVtxArray (TObjArray *InVtxArray)
 
virtual ~EdbShowAlg_GS ()
 
- Public Member Functions inherited from EdbShowAlg
void AddRecoShowerArray (EdbShowerP *shower)
 
 ClassDef (EdbShowAlg, 1)
 
Double_t DeltaR_NoPropagation (EdbSegP *s, EdbSegP *stest)
 
Double_t DeltaR_WithoutPropagation (EdbSegP *s, EdbSegP *stest)
 
Double_t DeltaR_WithPropagation (EdbSegP *s, EdbSegP *stest)
 
Double_t DeltaTheta (EdbSegP *s1, EdbSegP *s2)
 
Double_t DeltaThetaComponentwise (EdbSegP *s1, EdbSegP *s2)
 
Double_t DeltaThetaSingleAngles (EdbSegP *s1, EdbSegP *s2)
 
 EdbShowAlg ()
 
 EdbShowAlg (TString AlgName, Int_t AlgValue)
 
TString GetAlgName () const
 
Int_t GetAlgValue () const
 
Double_t GetMinimumDist (EdbSegP *seg1, EdbSegP *seg2)
 
TObjArray * GetRecoShowerArray () const
 
Int_t GetRecoShowerArrayN () const
 
EdbShowerPGetShower (Int_t i) const
 
Double_t GetSpatialDist (EdbSegP *s1, EdbSegP *s2)
 
void Help ()
 
Bool_t IsInConeTube (EdbSegP *sTest, EdbSegP *sStart, Double_t CylinderRadius, Double_t ConeAngle)
 
void Print ()
 
void PrintAll ()
 
void PrintMore ()
 
void PrintParameters ()
 
void PrintParametersShort ()
 
void SetActualAlgParameterset (Int_t ActualAlgParametersetNr)
 
void SetEdbPVRec (EdbPVRec *Ali)
 
void SetEdbPVRecPIDNumbers (Int_t FirstPlate_eAliPID, Int_t LastPlate_eAliPID, Int_t MiddlePlate_eAliPID, Int_t NumberPlate_eAliPID)
 
void SetInBTArray (TObjArray *InBTArray)
 
void SetParameter (Int_t parNr, Float_t par)
 
void SetParameters (Float_t *par)
 
void SetRecoShowerArray (TObjArray *RecoShowerArray)
 
void SetRecoShowerArrayN (Int_t RecoShowerArrayN)
 
void SetUseAliSub (Bool_t UseAliSub)
 
void Transform_eAli (EdbSegP *InitiatorBT, Float_t ExtractSize)
 
void UpdateShowerIDs ()
 
void UpdateShowerMetaData ()
 
virtual ~EdbShowAlg ()
 

Private Attributes

TObjArray * eInVtxArray
 
Int_t eInVtxArrayN
 
TH1F * h_GammaChi2
 

Additional Inherited Members

- Protected Member Functions inherited from EdbShowAlg
void Set0 ()
 
- Protected Attributes inherited from EdbShowAlg
Int_t eActualAlgParametersetNr
 
TString eAlgName
 
Int_t eAlgValue
 
EdbPVReceAli
 
EdbPVReceAli_Sub
 
Int_t eAli_SubNpat
 
Int_t eAliNpat
 
Int_t eFirstPlate_eAliPID
 
TObjArray * eInBTArray
 
Int_t eInBTArrayN
 
Int_t eLastPlate_eAliPID
 
Int_t eMiddlePlate_eAliPID
 
Int_t eNumberPlate_eAliPID
 
Int_t eParaN
 
TString eParaString [10]
 
Float_t eParaValue [10]
 
EdbShowerPeRecoShower
 
TObjArray * eRecoShowerArray
 
Int_t eRecoShowerArrayN
 
Int_t eUseAliSub
 

Constructor & Destructor Documentation

◆ EdbShowAlg_GS() [1/2]

EdbShowAlg_GS::EdbShowAlg_GS ( )

◆ EdbShowAlg_GS() [2/2]

EdbShowAlg_GS::EdbShowAlg_GS ( EdbPVRec ali)
46 {
47  // Default Constructor with EdbPVRec object
48  Log(2,"EdbShowAlg_GS::EdbShowAlg_GS(EdbPVRec* ali)","EdbShowAlg_GS:: Default Constructor with EdbPVRec object");
49  // Default Constructor
50  Log(2,"EdbShowAlg_GS::EdbShowAlg_GS","EdbShowAlg_GS:: Default Constructor");
51 
52  // Reset all:
53  Set0();
54 
55  eAlgName="GS";
56  eAlgValue=999; // see default.par_SHOWREC for labeling (labeling identical with ShowRec program)
57 
58  eInVtxArrayN=0;
60 
61  eInBTArrayN=0;
63 
64  // Init with values according to GS Alg:
65  Init();
66 
67  this->SetEdbPVRec(ali);
68 
69 }
bool Log(int level, const char *location, const char *fmt,...)
Definition: EdbLog.cxx:75
TObjArray * eInVtxArray
Definition: EdbShowAlg_GS.h:44
Int_t eInVtxArrayN
Definition: EdbShowAlg_GS.h:45
void Init()
Definition: EdbShowAlg_GS.cxx:80
TString eAlgName
Definition: EdbShowAlg.h:49
void Set0()
Definition: EdbShowAlg.cxx:53
Int_t eAlgValue
Definition: EdbShowAlg.h:50
void SetEdbPVRec(EdbPVRec *Ali)
Definition: EdbShowAlg.h:100
Int_t eInBTArrayN
Definition: EdbShowAlg.h:63
TObjArray * eInBTArray
Definition: EdbShowAlg.h:62
#define NULL
Definition: nidaqmx.h:84

◆ ~EdbShowAlg_GS()

EdbShowAlg_GS::~EdbShowAlg_GS ( )
virtual
73 {
74  // Default Destructor
75  Log(2,"EdbShowAlg_GS::~EdbShowAlg_GS","EdbShowAlg_GS:: Default Destructor");
76 }

Member Function Documentation

◆ CalcIP() [1/2]

double EdbShowAlg_GS::CalcIP ( EdbSegP s,
double  x,
double  y,
double  z 
)
688  {
689  // Calculate IP between a given segment and a given x,y,z.
690  // return the IP value.
691  double ax = s->TX();
692  double ay = s->TY();
693  double bx = s->X()-ax*s->Z();
694  double by = s->Y()-ay*s->Z();
695  double a;
696  double r;
697  double xx,yy,zz;
698  a = (ax*(x-bx)+ay*(y-by)+1.*(z-0.))/(ax*ax+ay*ay+1.);
699  xx = bx +ax*a;
700  yy = by +ay*a;
701  zz = 0. +1.*a;
702  r = sqrt((xx-x)*(xx-x)+(yy-y)*(yy-y)+(zz-z)*(zz-z));
703  return r;
704 }
void a()
Definition: check_aligned.C:59
Float_t TX() const
Definition: EdbSegP.h:172
Float_t X() const
Definition: EdbSegP.h:170
Float_t Z() const
Definition: EdbSegP.h:150
Float_t Y() const
Definition: EdbSegP.h:171
Float_t TY() const
Definition: EdbSegP.h:173
void r(int rid=2)
Definition: test.C:201
EdbSegP * s
Definition: tlg2pattern.C:32

◆ CalcIP() [2/2]

double EdbShowAlg_GS::CalcIP ( EdbSegP s,
EdbVertex v 
)
708  {
709  // calculate IP between a given segment and a given vertex.
710  // return the IP value.
711  // this is used for IP cut.
712 
713  // if vertex is not given, use the selected vertex.
714  // if(v==NULL) v=gEDA->GetSelectedVertex();
715 
716  if (v==NULL) return -1.;
717  return CalcIP(s, v->X(), v->Y(), v->Z());
718 }
Double_t CalcIP(EdbSegP *s, double x, double y, double z)
Definition: EdbShowAlg_GS.cxx:688
Float_t X() const
Definition: EdbVertex.h:130
Float_t Z() const
Definition: EdbVertex.h:132
Float_t Y() const
Definition: EdbVertex.h:131

◆ CheckCleanPairs()

TObjArray * EdbShowAlg_GS::CheckCleanPairs ( EdbSegP InBT,
TObjArray *  RecoShowerArrayFromFindPairs 
)
614 {
615 
616  cout << "EdbShowAlg_GS::CheckCleanPairs" <<endl;
617  if (NULL==RecoShowerArray) return NULL;
618 
619  TObjArray* NewRecoShowerArray= new TObjArray(999);
620  int NewRecoShowerArrayN=0;
621 
622  EdbShowerP* TrackPair1=NULL;
623  EdbShowerP* TrackPair2=NULL;
624 
625  int ntrack=RecoShowerArray->GetEntriesFast();
626  cout << ntrack << endl;
627 
628  for (Int_t pat_one_cnt=0; pat_one_cnt<ntrack; ++pat_one_cnt) {
629  // (if -1) then the last one is not checked
630  TrackPair1=(EdbShowerP*)RecoShowerArray->At(pat_one_cnt);
631  bool taketrack1=true;
632  for (Int_t pat_two_cnt=pat_one_cnt; pat_two_cnt<ntrack; ++pat_two_cnt) {
633  //if only one track at all take it anyway:
634  if (ntrack==1) continue;
635  TrackPair2=(EdbShowerP*)RecoShowerArray->At(pat_two_cnt);
636  // cout << pat_one_cnt << " " << pat_two_cnt << endl;
637  // cout << pat_one_cnt << " " << pat_two_cnt << endl;
638 
639  // Check if track1 has a track before (smaller Z value)
640  // and if so, if dtheta is smaller than 0.1:
641  // If both is so, then track 1 is NOT taken for final array:
642  // cout << TrackPair1->Z() << " " << TrackPair2->Z() << endl;
643  EdbSegP* s1=(EdbSegP* )TrackPair1->GetSegment(0);
644  // s1->PrintNice();
645  EdbSegP* s2=(EdbSegP* )TrackPair2->GetSegment(0);
646  // s2->PrintNice();
647  if (TrackPair1->Z()>TrackPair2->Z() && DeltaThetaSingleAngles((EdbSegP*)s1,(EdbSegP*)s2)<0.05) taketrack1=false;
648  if (!taketrack1) break;
649  }
650 
651  if (!taketrack1) continue;
652 
653  // Add TrackPair1
654  NewRecoShowerArray->Add(TrackPair1);
655  ++NewRecoShowerArrayN;
656 
657 
658  // Print
659  if (TrackPair1->N()<2) cout << "This Track has only ONE entry" << endl;
660 
661  EdbSegP* s1=(EdbSegP* )TrackPair1->GetSegment(0);
662  EdbSegP* s2=(EdbSegP* )TrackPair1->GetSegment(1);
663  if (gEDBDEBUGLEVEL>2) {
664  s1->PrintNice();
665  s2->PrintNice();
666  cout << s1->MCEvt() << " " << s2->MCEvt() << " " << s1->Flag() << " " << s2->Flag() << " " << s1->P() << " " << s2->P() << " " << " " << s1->Z() << " " << s2->Z() << " " <<endl;
667  cout<< "--------"<<endl;
668  }
669  }
670  cout << " From " <<ntrack << " originally, there are now after Zposition and overlap cuts: " << NewRecoShowerArrayN << " left." << endl;
671 
672  return NewRecoShowerArray;
673 }
TObjArray * RecoShowerArray
Definition: Shower_E_FromShowerRoot.C:12
Definition: EdbSegP.h:18
Float_t P() const
Definition: EdbSegP.h:149
void PrintNice() const
Definition: EdbSegP.cxx:418
Int_t MCEvt() const
Definition: EdbSegP.h:142
Int_t Flag() const
Definition: EdbSegP.h:146
Double_t DeltaThetaSingleAngles(EdbSegP *s1, EdbSegP *s2)
Definition: EdbShowAlg.cxx:455
Definition: EdbShowerP.h:28
Int_t N() const
Definition: EdbShowerP.h:412
EdbSegP * GetSegment(int i) const
Definition: EdbShowerP.h:435
gEDBDEBUGLEVEL
Definition: energy.C:7
EdbSegP * s1
Definition: tlg2pattern.C:30
EdbSegP * s2
Definition: tlg2pattern.C:31

◆ CheckPairDuplications()

Bool_t EdbShowAlg_GS::CheckPairDuplications ( Int_t  SegPID,
Int_t  SegID,
Int_t  Seg2PID,
Int_t  Seg2ID,
TArrayI *  SegmentPIDArray,
TArrayI *  SegmentIDArray,
TArrayI *  Segment2PIDArray,
TArrayI *  Segment2IDArray,
int  RecoShowerArrayN 
)
237 {
238  if (gEDBDEBUGLEVEL>2) cout << "EdbShowAlg_GS::CheckPairDuplications for Pair (" << SegPID << "," << SegID<< ";"<< Seg2PID << "," << Seg2ID << ") in RecoShowerArrayN=" << RecoShowerArrayN << endl;
239 
240  for (Int_t i=0; i<RecoShowerArrayN; i++) {
241  // PID and ID of Seg and Seg2 to be exchanged for duplications
242  if (gEDBDEBUGLEVEL>3) cout << "EdbShowAlg_GS::CheckPairDuplications in eRecoShowerArray: i="<< i<<" (" << Segment2PIDArray->At(i) << "," << Segment2IDArray->At(i)<< ";"<< SegmentPIDArray->At(i) << "," << SegmentIDArray->At(i) << "). "<< endl;
243  if ( SegPID==Segment2PIDArray->At(i) && Seg2PID==SegmentPIDArray->At(i) && SegID==Segment2IDArray->At(i) && Seg2ID==SegmentIDArray->At(i)) {
244  if (gEDBDEBUGLEVEL>3) cout << "Found duplication for i="<< i<< ".. return true"<<endl;
245  return kTRUE;
246  }
247  }
248  if (gEDBDEBUGLEVEL>2) cout << "EdbShowAlg_GS::CheckPairDuplications ... no duplication found."<<endl;
249  return kFALSE;
250 }

◆ ClassDef()

EdbShowAlg_GS::ClassDef ( EdbShowAlg_GS  ,
 
)

◆ Convert_InVtxArray_To_InBTArray()

void EdbShowAlg_GS::Convert_InVtxArray_To_InBTArray ( )
150 {
151  Log(2,"EdbShowAlg_GS::EdbShowAlg_GS","Convert_InVtxArray_To_InBTArray()");
152  if (eInBTArray!=NULL) {
153  eInBTArray->Clear();
154  cout << " eInBTArray->Clear();" << endl;
155  }
156  if (eInBTArray==NULL) {
157  eInBTArray = new TObjArray();
158  cout << " eInBTArray = new TObjArray()" << endl;
159  }
160 
161  if (eInVtxArray==NULL) cout << "NO eInVtxArray " << endl;
162  EdbVertex* vtx;
163  cout << "EdbVertex* vtx;" << endl;
164  cout << eInVtxArrayN << endl;
165  cout << "now eInVtxArray->Print();" << endl;
166  eInVtxArray->Print();
167  cout << eInVtxArray->GetEntries() << endl;
168 
169  Log(2,"EdbShowAlg_GS::EdbShowAlg_GS","Convert_InVtxArray_To_InBTArray()...start loop");
170  for (Int_t i=0; i<eInVtxArrayN; i++) {
171  vtx= (EdbVertex*)eInVtxArray->At(i);
172  //vtx->Print();
173  EdbSegP* seg = new EdbSegP(i,vtx->X(),vtx->Y(),0,0);
174  seg->SetZ(vtx->Z());
175  // vtx can have by default initialization MCEvt==0,
176  // this we dont want, in case: we set MCEvt of vtx from 0 to -999
177  if (vtx->MCEvt()==0) vtx->SetMC(-999);
178  seg->SetMC(vtx->MCEvt(),vtx->MCEvt());
179  seg->SetFlag(0);
180  //seg->Print();
181  eInBTArray->Add(seg);
182  }
183 
184  eInBTArrayN=eInBTArray->GetEntries();
185  if (gEDBDEBUGLEVEL>2) cout << "EdbShowAlg_GS::Convert_InVtxArray_To_InBTArray Converted " << eInBTArrayN << "InVtx to InBT." << endl;
186  Log(2,"EdbShowAlg_GS::EdbShowAlg_GS","Convert_InVtxArray_To_InBTArray()...done.");
187  return;
188 }
void SetZ(float z)
Definition: EdbSegP.h:122
void SetMC(int mEvt, int mTrack)
Definition: EdbSegP.h:138
void SetFlag(int flag)
Definition: EdbSegP.h:127
Definition: EdbVertex.h:68
void SetMC(int mEvt=0)
Definition: EdbVertex.h:159
Int_t MCEvt() const
Definition: EdbVertex.h:125

◆ Convert_TracksArray_To_InBTArray()

void EdbShowAlg_GS::Convert_TracksArray_To_InBTArray ( )
194 {
195  Log(2,"EdbShowAlg_GS::EdbShowAlg_GS","Convert_TracksArray_To_InBTArray()");
196  if (eInBTArray!=NULL) {
197  eInBTArray->Clear();
198  cout << " eInBTArray->Clear();" << endl;
199  }
200  if (eInBTArray==NULL) {
201  eInBTArray = new TObjArray();
202  cout << " eInBTArray = new TObjArray()" << endl;
203  }
204 
205  TObjArray* tracks=eAli->eTracks;
206  Int_t ntracks=tracks->GetEntries();
207  EdbSegP* sg=0;
208  EdbTrackP* track=0;
209 
210  Log(2,"EdbShowAlg_GS::EdbShowAlg_GS","Convert_TracksArray_To_InBTArray()...start loop");
211  for (Int_t i=0; i<ntracks; i++) {
212  track= (EdbTrackP*)tracks->At(i);
213  sg= (EdbSegP*)track->GetSegmentFirst();
214  //vtx->Print();
215  EdbSegP* seg = new EdbSegP(i,sg->X(),sg->Y(),0,0);
216  seg->SetZ(sg->Z());
217  seg->SetTX(sg->TX());
218  seg->SetTY(sg->TY());
219  seg->SetMC(sg->MCEvt(),sg->MCEvt());
220  seg->SetFlag(0);
221  seg->PropagateTo(seg->Z()-1300); // propagate one plate downstream;
222  //seg->Print();
223  eInBTArray->Add(seg);
224  }
225 
226  eInBTArrayN=eInBTArray->GetEntries();
227  if (gEDBDEBUGLEVEL>1) cout << "EdbShowAlg_GS::Convert_TracksArray_To_InBTArray Converted " << eInBTArrayN << " to InBT." << endl;
228  Log(2,"EdbShowAlg_GS::EdbShowAlg_GS","Convert_TracksArray_To_InBTArray()...done.");
229  return;
230 }
TObjArray * eTracks
Definition: EdbPVRec.h:161
void PropagateTo(float z)
Definition: EdbSegP.cxx:292
void SetTX(Float_t tx)
Definition: EdbSegP.h:176
void SetTY(Float_t ty)
Definition: EdbSegP.h:177
EdbPVRec * eAli
Definition: EdbShowAlg.h:59
Definition: EdbPattern.h:118
Definition: bitview.h:14
TTree * tracks
Definition: check_tr.C:19

◆ Execute()

void EdbShowAlg_GS::Execute ( )
virtual

Reimplemented from EdbShowAlg.

256 {
257  Log(2,"EdbShowAlg_GS::Execute","Execute()");
258  cout << "EdbShowAlg_GS::Execute()" << endl;
259  Log(2,"EdbShowAlg_GS::Execute","Execute() DOING MAIN SHOWER RECONSTRUCTION HERE");
260 
261  if (!eAli) cout << "no eAli pointer!" << endl;
262 
263  if (!eInBTArray) {
264  Log(2,"EdbShowAlg_GS::Execute","Execute() No eInBTArray.... Check for eInVtxArray:");
265  Log(2,"EdbShowAlg_GS::Execute","Execute() Check if there is an Convert_InVtxArray_To_InBTArray():");
266  // Check if there is an Convert_InVtxArray_To_InBTArray();
267  if (!eInVtxArray) {
268  Log(2,"EdbShowAlg_GS::Execute","Execute() No eInVtxArray. RETURN.");
269 
270  cout << "Check if EdbPVrecObject has Tracks..." << endl;
271  if (NULL==eAli->eTracks) cout << " EdbPVrecObject has no Tracks" << endl;
272  cout << "... yes, seems so, so take the track starting segments, propagate their first segemnts one plate downstream" << endl;
273  cout << "and use this as invtx array...." << endl;
275  return;
276  }
277  cout << "eInVtxArray there. Converting to InBTArray() now. "<< endl;
279  }
280 
281  EdbSegP* InBT=NULL;
282  Int_t STEP=-1;
284  if (gEDBDEBUGLEVEL>1) cout << "EdbShowAlg_GS::Execute STEP for patternloop direction = " << STEP << endl;
285  if (gEDBDEBUGLEVEL>1) cout << "EdbShowAlg_GS::Execute eRecoShowerArrayN = " << eRecoShowerArrayN << endl;
286 
287  //--- Loop over InBTs:
288  if (gEDBDEBUGLEVEL>1) cout << "EdbShowAlg_GS::Execute Loop over InBTs:" << endl;
289 
290  // Since eInBTArray is filled in ascending ordering by zpositon
291  // We use the descending loop to begin with BT with lowest z first.
292  // InBT doest have to be necessary a real BaseTrack, it can also be a vertex (with its positions and angle zero) !!!
293  for (Int_t i=eInBTArrayN-1; i>=0; --i) {
294 
295  // CounterOutPut
296  if (gEDBDEBUGLEVEL==2) if ((i%1)==0) ;
297  // cout << eInBTArrayN <<" InBT in total, still to do:"<<Form("%4d",i)<< "\r\r\r\r"<<flush;
298  cout << "EdbShowAlg_GS::Execute ----------------------------------------------"<< endl;
299  cout << "EdbShowAlg_GS::Execute InBT in total, still to do: " << i << endl;
300 
301  // Get InitiatorBT from eInBTArray InBT
302  InBT=(EdbSegP*)eInBTArray->At(i);
303  cout << "EdbShowAlg_GS::Execute InBT->XYZ() " << InBT->X() << " " << InBT->Y() << " " << InBT->Z() << " " << endl;
304 
305  //-----------------------------------
306  // 1) Make local_gAli with cut parameters:
307  //-----------------------------------
308  // Transform (make size smaller, extract only events having same MC) the eAli object:
309  Transform_eAli(InBT,999999);
310  if (gEDBDEBUGLEVEL>2) eAli_Sub->Print();
311 
312  //-----------------------------------
313  // 2) Find Pairs
314  //-----------------------------------
315  TObjArray* Pairs = FindPairs(InBT,eAli_Sub);
316  cout << "EdbShowAlg_GS::Execute TObjArray* Pairs = FindPairs(InBT,eAli_Sub); done." << endl;
317  cout << "EdbShowAlg_GS::Execute Found Pair Entries= " << Pairs->GetEntries() << ". Do clean of pairings now." << endl;
318 
319 
320  //-----------------------------------
321  // 2) Add Clean Pairs to eRecoShowerArray
322  //-----------------------------------
323  TObjArray* CleanPairs = CheckCleanPairs( InBT, Pairs);
324  if (gEDBDEBUGLEVEL>2) cout << "EdbShowAlg_GS::Execute CleanPairs->Print(): "<< endl;
325  if (gEDBDEBUGLEVEL>2) CleanPairs->Print();
326  for (int j=0; j<CleanPairs->GetEntries(); ++j) {
327  //cout << "j= " << j << endl;
328  EdbShowerP* pair=(EdbShowerP*)CleanPairs->At(j);
329  eRecoShowerArray->Add(pair);
331  }
332 
333  cout << "EdbShowAlg_GS::Execute eRecoShowerArrayN now: " << eRecoShowerArrayN << endl;
334  cout << "EdbShowAlg_GS::Execute ----------------------------------------------"<< endl;
335  } // Loop over InBT
336 
337 
338  if (gEDBDEBUGLEVEL>2) {
339  cout << "EdbShowAlg_GS::Execute Print Loop for N shower pairsegments: = " << eRecoShowerArrayN << endl;
340  for (Int_t i=0; i <eRecoShowerArray->GetEntries(); i++) {
341  EdbShowerP* sh = (EdbShowerP* )eRecoShowerArray->At(i);
342  sh->Print();
343  sh->PrintSegments();
344  }
345  cout << "EdbShowAlg_GS::Execute Print Loop for N shower pairsegments ..... done."<< endl;
346  }
347 
348  cout << "EdbShowAlg_GS::Execute Found eRecoShowerArrayN Gamma Pairing Candidates: = " << eRecoShowerArrayN << endl;
349 
350  Log(2,"EdbShowAlg_GS::Execute","Execute()...done.");
351  return;
352 }
void Print() const
Definition: EdbPattern.cxx:1853
TObjArray * CheckCleanPairs(EdbSegP *InBT, TObjArray *RecoShowerArrayFromFindPairs)
Definition: EdbShowAlg_GS.cxx:613
void Convert_TracksArray_To_InBTArray()
Definition: EdbShowAlg_GS.cxx:193
TObjArray * FindPairs(EdbSegP *InBT, EdbPVRec *eAli_Sub)
Definition: EdbShowAlg_GS.cxx:359
void Convert_InVtxArray_To_InBTArray()
Definition: EdbShowAlg_GS.cxx:149
void Transform_eAli(EdbSegP *InitiatorBT, Float_t ExtractSize)
Definition: EdbShowAlg.cxx:107
TObjArray * eRecoShowerArray
Definition: EdbShowAlg.h:70
Int_t eFirstPlate_eAliPID
Definition: EdbShowAlg.h:65
Int_t eLastPlate_eAliPID
Definition: EdbShowAlg.h:66
Int_t eRecoShowerArrayN
Definition: EdbShowAlg.h:71
EdbPVRec * eAli_Sub
Definition: EdbShowAlg.h:74
void Print()
-— EXPERIMENTAL STATUS —
Definition: EdbShowerP.cxx:2333
void PrintSegments()
Definition: EdbShowerP.cxx:2367

◆ Finalize()

void EdbShowAlg_GS::Finalize ( )
virtual

Reimplemented from EdbShowAlg.

682 {
683  // do nothing yet.
684 }

◆ FindPairs()

TObjArray * EdbShowAlg_GS::FindPairs ( EdbSegP InBT,
EdbPVRec eAli_Sub 
)
360 {
361  Log(3,"EdbShowAlg_GS::FindPairs","FindPairs()");
362 
363  TObjArray* RecoShowerArray= new TObjArray(99);
364  Int_t RecoShowerArrayN=0;
365  TArrayI* SegmentPIDArray = new TArrayI(9999);
366  TArrayI* SegmentIDArray = new TArrayI(9999);
367  TArrayI* Segment2PIDArray = new TArrayI(9999);
368  TArrayI* Segment2IDArray = new TArrayI(9999);
369  EdbSegP* Segment=NULL;
370  EdbSegP* Segment2=NULL;
371  Float_t x_av,y_av,z_av,tx_av,ty_av,distZ;
372  EdbSegP* Segment_Sum=new EdbSegP(0,0,0,0,0,0);
373 
374  Float_t IP_Pair_To_InBT=0;
375  Int_t IP_Pair_To_InBT_SegSmaller=0;
376 
377  EdbSegP* InBT=NULL;
378  if (NULL==InitiatorBT) {
379  InBT= new EdbSegP();
381  InBT->SetX(pat->X());
382  InBT->SetY(pat->Y());
383  InBT->SetZ(pat->Z());
384  InBT->SetTX(0);
385  InBT->SetTY(0);
386  InBT->SetMC(-999,-999);
387  cout << "WARNING EdbShowAlg_GS::FindPairs InBT==NULL. Create a dummy InBT:" << endl;
388  InBT->PrintNice();
389  }
390  else {
391  InBT=InitiatorBT;
392  }
393 
394 
395 // from ShowRec.cpp:
396 // CUT_PARAMETER[0]=cut_gs_cut_dip;
397 // CUT_PARAMETER[1]=cut_gs_cut_dmin;
398 // CUT_PARAMETER[2]=cut_gs_cut_dr;
399 // CUT_PARAMETER[3]=cut_gs_cut_dz;
400 // CUT_PARAMETER[4]=cut_gs_cut_dtheta;
401 // CUT_PARAMETER[5]=cut_gs_cut_piddiff;
402 // CUT_PARAMETER[6]=cut_gs_cut_oppositeflag;
403 
404 
405  //-----------------------------------
406  // 2) Loop over (whole) eAli, check BT for Cuts
407  // eAli_Sub
408  // Loop structure:
409  // Loop over plates [0..NPatterns-1] for pattern one
410  // Loop over plates [0..NPatterns-1] for pattern two
411  // Take only plate pairings for |PID diff|<=3
412  // Loop over all BT of pattern one
413  // Loop over all BT of pattern two
414 
415  // This doesnt yet distinguish the FIRSTPLATE, MIDDLEPLATE, LATPLATE,NUMBERPLATES
416  // labelling, this will be built in later....
417  //-----------------------------------
418 
419  Int_t npat=eAli_Sub->Npatterns();
420  Int_t pat_one_bt_cnt_max,pat_two_bt_cnt_max=0;
421  EdbPattern* pat_one=0;
422  EdbPattern* pat_two=0;
423 
424  for (Int_t pat_one_cnt=0; pat_one_cnt<npat; ++pat_one_cnt) {
425  pat_one=(EdbPattern*)eAli_Sub->GetPattern(pat_one_cnt);
426  pat_one_bt_cnt_max=eAli_Sub->GetPattern(pat_one_cnt)->GetN();
427 
428  // Check if dist Z to vtx (BT) is ok:
429  distZ=pat_one->Z()-InBT->Z();
430  if (gEDBDEBUGLEVEL>4) {
431  cout << "EdbShowerAlg_GS::FindPairs Check if dist Z to vtx (BT) is ok: distZ=" << distZ << endl;
432  }
433  if (distZ>eParaValue[3]) continue;
434  if (distZ<0) continue;
435 
436  if (gEDBDEBUGLEVEL>2) cout << "Searching patterns: pat_one_cnt=" << pat_one_cnt << " pat_one->Z() = " << pat_one->Z() << " pat_one_bt_cnt_max= "<< pat_one_bt_cnt_max <<endl;
437 
438 
439  for (Int_t pat_two_cnt=0; pat_two_cnt<npat; ++pat_two_cnt) {
440 
441  // Now apply cut conditions: GS GAMMA SEARCH Alg:
442  if (TMath::Abs(pat_one_cnt-pat_two_cnt)>eParaValue[5]) continue;
443 
444  pat_two=(EdbPattern*)eAli_Sub->GetPattern(pat_two_cnt);
445  pat_two_bt_cnt_max=eAli_Sub->GetPattern(pat_two_cnt)->GetN();
446 
447 
448  if (gEDBDEBUGLEVEL>2) cout << " Searching patterns: pat_two_cnt=" << pat_two_cnt << " pat_two->Z() = " << pat_two->Z() << " pat_two_bt_cnt_max= "<< pat_two_bt_cnt_max <<endl;
449 
450  for (Int_t pat_one_bt_cnt=0; pat_one_bt_cnt<pat_one_bt_cnt_max; ++pat_one_bt_cnt) {
451  Segment = (EdbSegP*)pat_one->GetSegment(pat_one_bt_cnt);
452 
453  Float_t IP_Pair_To_InBT_Seg =CalcIP(Segment, InBT->X(),InBT->Y(),InBT->Z());
454  // Check if IP to vtx (BT) is ok:
455  if (IP_Pair_To_InBT>eParaValue[0]) continue;
456 
457 
458  for (Int_t pat_two_bt_cnt=0; pat_two_bt_cnt<pat_two_bt_cnt_max; ++pat_two_bt_cnt) {
459 
460  Segment2 = (EdbSegP*)pat_two->GetSegment(pat_two_bt_cnt);
461  if (Segment2==Segment) continue;
462  if (Segment2->ID()==Segment->ID()&&Segment2->PID()==Segment->PID()) continue;
463 
464  // At first: Check for already duplicated pairings:
465  if (CheckPairDuplications(Segment->PID(),Segment->ID(),Segment2->PID(),Segment2->ID(), SegmentPIDArray,SegmentIDArray,Segment2PIDArray,Segment2IDArray, RecoShowerArrayN)) continue;
466 
467  // Now apply cut conditions: GS GAMMA SEARCH Alg --------------------
468  // if InBT is flagged as MC InBT, take care that only BG or same MC basetracks are taken:
469  if (InBT->MCEvt()>0) if (Segment->MCEvt()>0&&Segment2->MCEvt()>0) if (Segment->MCEvt()!=Segment2->MCEvt()) continue;
470  if (InBT->MCEvt()>0) if (Segment->MCEvt()>0&&Segment2->MCEvt()>0) if (Segment->MCEvt()!=InBT->MCEvt()) continue;
471  if (InBT->MCEvt()>0) if (Segment->MCEvt()>0&&Segment2->MCEvt()>0) if (Segment2->MCEvt()!=InBT->MCEvt()) continue;
472 
473  // In case of two MC events, check for e+ e- pairs
474  // Do this ONLY IF parameter eParaValue[6] is set to choose different Flag() pairs:
475  if (InBT->MCEvt()>0 && eParaValue[6]==1) {
476  if (Segment->MCEvt()>0&&Segment2->MCEvt()>0) {
477  if ((Segment2->Flag()+Segment->Flag())!=0) continue;
478  }
479  }
480 
481  // a) Check dR between tracks:
482  if (DeltaR_WithPropagation(Segment,Segment2)>eParaValue[2]) continue;
483  // b) Check dT between tracks:
484  if (DeltaThetaSingleAngles(Segment,Segment2)>eParaValue[4]) continue;
485  // c) Check dMinDist between tracks:
486  if (GetMinimumDist(Segment,Segment2)>eParaValue[1]) continue;
487 
488  // d) Check if dist Z to vtx (BT) is ok:
489  distZ=Segment->Z()-InBT->Z();
490  if (distZ>eParaValue[3]) continue;
491 
492  x_av=Segment2->X()+(Segment->X()-Segment2->X())/2.0;
493  y_av=Segment2->Y()+(Segment->Y()-Segment2->Y())/2.0;
494  z_av=Segment2->Z()+(Segment->Z()-Segment2->Z())/2.0;
495  tx_av=Segment2->TX()+(Segment->TX()-Segment2->TX())/2.0;
496  ty_av=Segment2->TY()+(Segment->TY()-Segment2->TY())/2.0;
497  Segment_Sum->SetX(x_av);
498  Segment_Sum->SetY(y_av);
499  Segment_Sum->SetTX(tx_av);
500  Segment_Sum->SetTY(ty_av);
501  Segment_Sum->SetZ(z_av);
502 
503 
504  // Check if IP to vtx (BT) is ok:
505  Float_t IP_Pair_To_InBT_Seg2 =CalcIP(Segment2, InBT->X(),InBT->Y(),InBT->Z());
506  if (IP_Pair_To_InBT_Seg2>eParaValue[0]) continue;
507 
508 
509 
510  Float_t IP_Pair_To_InBT_SegSum=CalcIP(Segment_Sum, InBT->X(),InBT->Y(),InBT->Z());
511  Float_t IP_Pair_To_InBT_SegSmaller=0;
512 
513 // cout << "IP_Pair_To_InBT_Seg = " << IP_Pair_To_InBT_Seg << endl;
514 // cout << "IP_Pair_To_InBT_Seg2 = " << IP_Pair_To_InBT_Seg2 << endl;
515 // cout << "IP_Pair_To_InBT_SegSum= " << IP_Pair_To_InBT_SegSum << endl;
516 
517  // Save the segment which has smaller IP, this will be the first BT in the RecoShower
518  if ( IP_Pair_To_InBT_Seg>IP_Pair_To_InBT_Seg2 ) {
519  IP_Pair_To_InBT_SegSmaller=0; // take Segment first
520  }
521  else {
522  IP_Pair_To_InBT_SegSmaller=1; // take Segment2 first
523  }
524 
525 
526 
527  // f) Check if this is not a possible fake doublet which is
528  // sometimes caused by view overlap in the scanning:
529  // in the EdbPVRQuality class this will be done at start for the whole
530  // PVR object so this will be later on obsolete.
531  if (IsPossibleFakeDoublet(Segment,Segment2) ) continue;
532  //
533  // end of cut conditions: GS GAMMA SEARCH Alg --------------------
534  //
535 
536 
537  // Chi2 Variable:
538  Float_t dminDist=GetMinimumDist(Segment,Segment2);
539  Float_t dtheta=DeltaThetaSingleAngles(Segment,Segment2);
540  // preliminary mean and sigma values (taken from 1Gev, with highpur cutset), to be checked if they are the same for 0.5,1,2,4 GeV Photons
541  Float_t GammaChi2 = ((IP_Pair_To_InBT-80)*(IP_Pair_To_InBT-80)/60/60)+((dminDist-3.5)*(dminDist-3.5)/4.7/4.7)+((dtheta-0.021)*(dtheta-0.021)/0.012/0.012);
542  GammaChi2=GammaChi2/3.0;
543 
544 
545  if (gEDBDEBUGLEVEL>2) {
546  cout << "EdbShowAlg_GS::FindPairs Pair (" << Segment->PID() << "," << Segment->ID()<< ";"<< Segment2->PID() << "," << Segment2->ID() << ") has passed all cuts w.r.t to InBT:" << endl;
547  cout << "EdbShowAlg_GS::FindPairs IP_Pair_To_InBT = " << IP_Pair_To_InBT << endl;
548  cout << "EdbShowAlg_GS::FindPairs (IP_Pair_To_InBT==0 can mean that there was no InVtxArray and Vtx has been extrapolated from tracks.)" << endl;
549  cout << "EdbShowAlg_GS::FindPairs GetMinimumDist(Segment,Segment2) = " << GetMinimumDist(Segment,Segment2) << endl;
550  cout << "EdbShowAlg_GS::FindPairs CalcIP(BetterSegment,InBT) = " << IP_Pair_To_InBT << endl;
551  cout << "EdbShowAlg_GS::FindPairs CalcIP(Segment_Sum,InBT) = " << IP_Pair_To_InBT_SegSum << endl;
552  cout << "EdbShowAlg_GS::FindPairs GetSpatialDist(Segment_Sum,InBT) = " << GetSpatialDist(Segment_Sum,InBT) << endl;
553  cout << "EdbShowAlg_GS::FindPairs GammaChi2 = " << GammaChi2 << endl;
554  }
555 
556 
557  // Add IDs and PIDs for storage
558  SegmentPIDArray->AddAt(Segment->PID(),RecoShowerArrayN);
559  SegmentIDArray->AddAt(Segment->ID(),RecoShowerArrayN);
560  Segment2PIDArray->AddAt(Segment2->PID(),RecoShowerArrayN);
561  Segment2IDArray->AddAt(Segment2->ID(),RecoShowerArrayN);
562 
563  // Create new EdbShowerP Object for storage;
564  EdbShowerP* RecoShower = new EdbShowerP();
565  if (IP_Pair_To_InBT_SegSmaller==0) {
566  RecoShower -> AddSegment(Segment);
567  RecoShower -> AddSegment(Segment2);
568  }
569  else {
570  RecoShower -> AddSegment(Segment2);
571  RecoShower -> AddSegment(Segment);
572  }
573  // Set X and Y and Z values: /Take lower Z of both BTs)
574  RecoShower->SetZ(TMath::Min(Segment->Z(),Segment2->Z()));
575  RecoShower->SetX(Segment_Sum->X());
576  RecoShower->SetY(Segment_Sum->Y());
577  RecoShower->SetTX(Segment_Sum->TX());
578  RecoShower->SetTY(Segment_Sum->TY());
579  RecoShower->SetMC(InBT->MCEvt(),InBT->MCEvt());
580  RecoShower->SetID(RecoShowerArrayN);
581  RecoShower->SetPID(Segment->PID());
582  //RecoShower ->PrintNice();
583 
584  if (gEDBDEBUGLEVEL>2) cout <<"------------"<< endl;
585 
586 
587  // Add Shower to interim Array:
588  RecoShowerArray->Add(RecoShower);
589  ++RecoShowerArrayN;
590 
591  } //for (Int_t pat_two_bt_cnt=0; ...
592  } //for (Int_t pat_one_bt_cnt=0; ...
593  } //for (Int_t pat_two_cnt=0; ...
594  } //for (Int_t pat_one_cnt=0; ...
595 
596 
597  cout << "EdbShowAlg_GS::FindPairs For the InBT/Vtx at __" << InBT << "__, we have found " << RecoShowerArray->GetEntries() << " compatible pairs in the PVRec volume." << endl;
598 
599  // Delete unnecessary objects:
600  // important, else memory overflow!!
601  delete SegmentPIDArray;
602  delete Segment2PIDArray;
603  delete SegmentIDArray;
604  delete Segment2IDArray;
605 
606  Log(3,"EdbShowAlg_GS::FindPairs","FindPairs()....done.");
607  return RecoShowerArray;
608 }
Int_t npat
Definition: Xi2HatStartScript.C:33
Definition: EdbPattern.h:280
EdbPattern * GetPatternZLowestHighest(Bool_t lowestZ=kTRUE) const
Definition: EdbPattern.cxx:1982
Int_t Npatterns() const
Definition: EdbPattern.h:380
EdbPattern * GetPattern(int id) const
Definition: EdbPattern.cxx:1887
void SetPID(int pid)
Definition: EdbSegP.h:126
void SetY(Float_t y)
Definition: EdbSegP.h:175
Int_t ID() const
Definition: EdbSegP.h:144
void SetID(int id)
Definition: EdbSegP.h:125
void SetX(Float_t x)
Definition: EdbSegP.h:174
Int_t PID() const
Definition: EdbSegP.h:145
Float_t Z() const
Definition: EdbPattern.h:87
Float_t Y() const
Definition: EdbPattern.h:86
Int_t GetN() const
Definition: EdbPattern.h:65
EdbSegP * GetSegment(int i) const
Definition: EdbPattern.h:66
Float_t X() const
Definition: EdbPattern.h:85
Bool_t IsPossibleFakeDoublet(EdbSegP *s1, EdbSegP *s2)
Definition: EdbShowAlg_GS.cxx:722
Bool_t CheckPairDuplications(Int_t SegPID, Int_t SegID, Int_t Seg2PID, Int_t Seg2ID, TArrayI *SegmentPIDArray, TArrayI *SegmentIDArray, TArrayI *Segment2PIDArray, TArrayI *Segment2IDArray, int RecoShowerArrayN)
Definition: EdbShowAlg_GS.cxx:236
Float_t eParaValue[10]
Definition: EdbShowAlg.h:52
Double_t GetSpatialDist(EdbSegP *s1, EdbSegP *s2)
Definition: EdbShowAlg.cxx:463
Double_t GetMinimumDist(EdbSegP *seg1, EdbSegP *seg2)
Definition: EdbShowAlg.cxx:1943
Double_t DeltaR_WithPropagation(EdbSegP *s, EdbSegP *stest)
Definition: EdbShowAlg.cxx:411

◆ GetInVtxArray()

TObjArray* EdbShowAlg_GS::GetInVtxArray ( ) const
inline
67  {
68  return eInVtxArray;
69  }

◆ GetInVtxArrayN()

Int_t EdbShowAlg_GS::GetInVtxArrayN ( ) const
inline
64  {
65  return eInVtxArrayN;
66  }

◆ Init()

void EdbShowAlg_GS::Init ( void  )
81 {
82  Log(2,"EdbShowAlg_GS::EdbShowAlg_GS","Init()");
83 
84  // Init with values according to GS Alg:
85  // Took over from "FindGamma.C" script I develoved before:
86  // IP CUT; this cut is used for the -better- IP of both BTs to Vertex/BT
87  eParaValue[0]=150;
88  eParaString[0]="PARA_GS_CUT_dIP";
89 
90  // min minDist.e between pair BTs
91  eParaValue[1]=40;
92  eParaString[1]="PARA_GS_CUT_dMin";
93 
94  // min dR between pair BTs
95  eParaValue[2]=60;
96  eParaString[2]="PARA_GS_CUT_dR";
97 
98  // max Z distance between pair BTs and Vertex/BT
99  eParaValue[3]=19000;
100  eParaString[3]="PARA_GS_CUT_dZ";
101 
102  // max Angle between pair BTs
103  eParaValue[4]=0.060;
104  eParaString[4]="PARA_GS_CUT_dtheta";
105 
106  // max plates difference between pair BTs
107  eParaValue[5]=1;
108  eParaString[5]="PARA_GS_CUT_PIDDIFF";
109 
110  // in MC case: check for opposite flag sign
111  eParaValue[6]=1;
112  eParaString[6]="PARA_GS_CUT_OPPOSITEFLAG";
113 
114  if (!eRecoShowerArray) eRecoShowerArray= new TObjArray();
116 
117  Log(2,"EdbShowAlg_GS::EdbShowAlg_GS","Init()...done.");
118  return;
119 }
TString eParaString[10]
Definition: EdbShowAlg.h:53

◆ Initialize()

void EdbShowAlg_GS::Initialize ( )
virtual

Reimplemented from EdbShowAlg.

125 {
126  Log(2,"EdbShowAlg_GS::EdbShowAlg_GS","Initialize()");
127  return;
128 }

◆ IsPossibleFakeDoublet()

Bool_t EdbShowAlg_GS::IsPossibleFakeDoublet ( EdbSegP s1,
EdbSegP s2 
)
722  {
723  if (TMath::Abs(s1->X()-s2->X())<1) return kTRUE; // minimum distance of 1micron
724  if (TMath::Abs(s1->Y()-s2->Y())<1) return kTRUE;// minimum distance of 1micron
725  if (TMath::Abs(s1->TX()-s2->TX())<0.005) return kTRUE;// minimum angle of 5mrad
726  if (TMath::Abs(s1->TY()-s2->TY())<0.005) return kTRUE;// minimum angle of 5mrad
727  return kFALSE;
728 }

◆ SetInVtx()

void EdbShowAlg_GS::SetInVtx ( EdbVertex vtx)
136 {
137  Log(2,"EdbShowAlg_GS::SetInVtx","SetInVtx()");
138  if (!eInVtxArray) eInVtxArray = new TObjArray();
139  eInVtxArray->Add(vtx);
140  ++eInVtxArrayN;
141  return;
142 }

◆ SetInVtxArray()

void EdbShowAlg_GS::SetInVtxArray ( TObjArray *  InVtxArray)
inline
57  {
58  eInVtxArray = InVtxArray;
59  eInVtxArrayN = eInVtxArray->GetEntries();
60  cout << eInVtxArrayN << " entries set"<<endl;
61  }

Member Data Documentation

◆ eInVtxArray

TObjArray* EdbShowAlg_GS::eInVtxArray
private

◆ eInVtxArrayN

Int_t EdbShowAlg_GS::eInVtxArrayN
private

◆ h_GammaChi2

TH1F* EdbShowAlg_GS::h_GammaChi2
private

The documentation for this class was generated from the following files: