FEDRA emulsion software from the OPERA Collaboration
Public Member Functions | Private Attributes | List of all members
EdbShowerAlg_GS Class Reference

#include <EdbShowerAlg.h>

Inheritance diagram for EdbShowerAlg_GS:
Inheritance graph
[legend]
Collaboration diagram for EdbShowerAlg_GS:
Collaboration graph
[legend]

Public Member Functions

void AddInVtx (EdbVertex *vtx)
 
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 CheckInput ()
 
Bool_t CheckPairDuplications (Int_t SegPID, Int_t SegID, Int_t Seg2PID, Int_t Seg2ID, TArrayI *SegmentPIDArray, TArrayI *SegmentIDArray, TArrayI *Segment2PIDArray, TArrayI *Segment2IDArray, Int_t RecoShowerArrayN)
 
 ClassDef (EdbShowerAlg_GS, 1)
 
void Convert_InVtxArray_To_InBTArray ()
 
void CreateANNPair ()
 
void CreateANNPlots ()
 
 EdbShowerAlg_GS ()
 
void Execute ()
 
void FillANNTree (TObjArray *RecoShowerArray, EdbSegP *Inbt)
 
void Finalize ()
 
TObjArray * FindPairs (EdbSegP *InBT, EdbPVRec *eAli_Sub)
 
TObjArray * FindPairsPreselected (EdbPVRec *eAli_Sub)
 
Bool_t FindPairsPreselected (EdbSegP *InitiatorBT, EdbPVRec *eAli_Sub)
 
TObjArray * GetInVtxArray () const
 
Int_t GetInVtxArrayN () const
 
void Init ()
 
void Initialize ()
 
Bool_t IsPossibleFakeDoublet (EdbSegP *s1, EdbSegP *s2)
 
void ReloadANNs (Int_t RecoMode)
 
TObjArray * SelectHighestPInMCArray (TObjArray *BTArray)
 
void Set0 ()
 
void SetCleanPairs (Bool_t CleanPairs)
 
void SetCutModeFull (Bool_t CutModeFull)
 
void SetFindPairsPreselected (Bool_t FindPairsPreselected)
 
void SetInVtx (EdbVertex *vtx)
 
void SetInVtxArray (TObjArray *InVtxArray)
 
void SetRecoMode (Int_t RecoMode)
 
virtual ~EdbShowerAlg_GS ()
 
- Public Member Functions inherited from EdbShowerAlg
void AddRecoShowerArray (TObjArray *RecoShowerArray)
 
EdbVertex * CalcVertex (TObjArray *segments)
 
 ClassDef (EdbShowerAlg, 1)
 
void Convert_EdbPVRec_To_InBTArray ()
 
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)
 
 EdbShowerAlg ()
 
 EdbShowerAlg (TString AlgName, Int_t AlgValue)
 
virtual void Execute ()
 
virtual void Finalize ()
 
TString GetAlgName () const
 
Int_t GetAlgValue () const
 
Double_t GetMinimumDist (EdbSegP *seg1, EdbSegP *seg2)
 
TObjArray * GetRecoShowerArray () const
 
Int_t GetRecoShowerArrayN () const
 
Double_t GetSpatialDist (EdbSegP *s1, EdbSegP *s2)
 
virtual void Initialize ()
 
Bool_t IsInConeTube (EdbSegP *sTest, EdbSegP *sStart, Double_t CylinderRadius, Double_t ConeAngle)
 
void Print ()
 
void PrintAll ()
 
void PrintMore ()
 
void PrintParameters ()
 
void PrintParametersShort ()
 
void PrintRecoShowerArray ()
 
void SetActualAlgParameterset (Int_t ActualAlgParametersetNr)
 
void SetDebug ()
 
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 parvalue)
 
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)
 
virtual ~EdbShowerAlg ()
 virtual constructor due to inherited class More...
 

Private Attributes

TMultiLayerPerceptron * eANNPair
 
TMultiLayerPerceptron * eANNPairCaseA
 
TMultiLayerPerceptron * eANNPairCaseB
 
TMultiLayerPerceptron * eANNPairCaseC
 
Float_t eANNPairCut [3]
 
TTree * eANNPairTree
 
Bool_t eCutModeFull
 
Bool_t eFindPairsPreselected
 
TObjArray * eInVtxArray
 
Int_t eInVtxArrayN
 
Bool_t eInVtxArraySet
 
Int_t eRecoMode
 
TArrayI * eSegment2IDArray
 
TArrayI * eSegment2PIDArray
 
TArrayI * eSegmentIDArray
 
TArrayI * eSegmentPIDArray
 
Bool_t eSetCleanPairs
 
Float_t eValueGSNN_var00
 
Float_t eValueGSNN_var01
 
Float_t eValueGSNN_var02
 
Float_t eValueGSNN_var03
 
Float_t eValueGSNN_var04
 
Float_t eValueGSNN_var05
 
Float_t eValueGSNN_var06
 
Float_t eValueGSNN_varInput
 
Float_t eValueGSNN_varOutput
 

Additional Inherited Members

- Protected Member Functions inherited from EdbShowerAlg
void Set0 ()
 
- Protected Attributes inherited from EdbShowerAlg
Int_t eActualAlgParametersetNr
 Used when more sets of same algorithm. More...
 
TString eAlgName
 
Int_t eAlgValue
 
EdbPVRec * eAli
 
EdbPVRec * eAli_Sub
 
Int_t eAli_SubNpat
 
Int_t eAliNpat
 
Bool_t eDebug
 
Int_t eFirstPlate_eAliPID
 
TObjArray * eInBTArray
 
Int_t eInBTArrayN
 
Int_t eLastPlate_eAliPID
 
Int_t eMiddlePlate_eAliPID
 
Int_t eNumberPlate_eAliPID
 
TString eParaString [10]
 
Float_t eParaValue [10]
 
EdbTrackP * eRecoShower
 
TObjArray * eRecoShowerArray
 
Int_t eRecoShowerArrayN
 
Int_t eUseAliSub
 

Constructor & Destructor Documentation

◆ EdbShowerAlg_GS()

EdbShowerAlg_GS::EdbShowerAlg_GS ( )

Default Constructor
This Shower Algorithm finds all compatible basetrack pairings that could
come from an ePlus-eMinus pair.
It needs an EdbPVRec* object.
Additionally -but not necessary- an array (or one) vertex, w.r.t. which the
basetracks should be looked at.
If vertex array is not given, then the basetracks of the volume will be
used to create interim helper vertices.
Usage:
Alg = new EdbShowerAlg_GS();
Alg->SetEdbPVRec(EdbPVRec* pvr);
Alg->SetInVtx(EdbVertex* vtx); // additionally...
Alg->Execute();
Alg->GetRecoShowerArray();
Returns the compatible Pair Segments (stored as EdbTrackP*).

703{
721
722 Log(2,"EdbShowerAlg_GS::EdbShowerAlg_GS","EdbShowerAlg_GS:: Default Constructor");
723
724 cout << "Howto do: " << endl;
725 cout << " Alg = new EdbShowerAlg_GS();" << endl;
726 cout << " Alg->SetEdbPVRec(EdbPVRec* pvr);" << endl;
727 cout << " Alg->SetInVtx(EdbVertex* vtx) // additionally... " << endl;
728 cout << " Alg->Execute();" << endl;
729 cout << " Alg->GetRecoShowerArray(); // .... Returns you the compatible Pair Segments (stored as EdbTrackP*)" << endl;
730 cout << " " << endl;
731
732
733 // Init with values according to GS Alg:
734 // Create essential objects:
735 Init();
736
737 // Reset all internal variables.
738 Set0();
739
740 // Create Neural Network for better fake BT rejection
741 CreateANNPair();
742}
Log
bool Log(int level, const char *location, const char *fmt,...)
Definition: EdbLog.cxx:75
EdbShowerAlg_GS::CreateANNPair
void CreateANNPair()
Definition: EdbShowerAlg.cxx:902
EdbShowerAlg_GS::Init
void Init()
Definition: EdbShowerAlg.cxx:811
EdbShowerAlg_GS::Set0
void Set0()
Definition: EdbShowerAlg.cxx:756

◆ ~EdbShowerAlg_GS()

EdbShowerAlg_GS::~EdbShowerAlg_GS ( )
virtual

Default Destructor

747{
749 Log(2,"EdbShowerAlg_GS::~EdbShowerAlg_GS","EdbShowerAlg_GS:: Default Destructor");
750}

Member Function Documentation

◆ AddInVtx()

void EdbShowerAlg_GS::AddInVtx ( EdbVertex *  vtx)
1051{
1052 Log(3,"EdbShowerAlg_GS::AddInVtx","AddInVtx()");
1053 eInVtxArray->Add(vtx);
1054 ++eInVtxArrayN;
1055 eInVtxArraySet=kTRUE;
1056 if (gEDBDEBUGLEVEL>2) {
1057 cout << "EdbShowerAlg_GS::AddInVtx(): Added vtx (XYZ,MC): " << vtx->X() << " " <<vtx->Y() << " " << vtx->Z() << " " << vtx->MCEvt() << " " << endl;
1058 cout << "EdbShowerAlg_GS::AddInVtx(): Added One Vertex. Now there are " << eInVtxArrayN << " InVtx stored." << endl;
1059 }
1060 Log(3,"EdbShowerAlg_GS::AddInVtx","AddInVtx()...done");
1061 return;
1062}
EdbShowerAlg_GS::eInVtxArraySet
Bool_t eInVtxArraySet
Definition: EdbShowerAlg.h:197
EdbShowerAlg_GS::eInVtxArrayN
Int_t eInVtxArrayN
Definition: EdbShowerAlg.h:196
EdbShowerAlg_GS::eInVtxArray
TObjArray * eInVtxArray
Definition: EdbShowerAlg.h:195
EdbVertex::X
Float_t X() const
Definition: EdbVertex.h:130
EdbVertex::MCEvt
Int_t MCEvt() const
Definition: EdbVertex.h:125
EdbVertex::Z
Float_t Z() const
Definition: EdbVertex.h:132
EdbVertex::Y
Float_t Y() const
Definition: EdbVertex.h:131
gEDBDEBUGLEVEL
gEDBDEBUGLEVEL
Definition: energy.C:7

◆ CalcIP() [1/2]

double EdbShowerAlg_GS::CalcIP ( EdbSegP *  s,
double  x,
double  y,
double  z 
)

Calculate IP between a given segment and a given x,y,z.
return the IP value.

2387 {
2390 double ax = s->TX();
2391 double ay = s->TY();
2392 double bx = s->X()-ax*s->Z();
2393 double by = s->Y()-ay*s->Z();
2394 double a;
2395 double r;
2396 double xx,yy,zz;
2397 a = (ax*(x-bx)+ay*(y-by)+1.*(z-0.))/(ax*ax+ay*ay+1.);
2398 xx = bx +ax*a;
2399 yy = by +ay*a;
2400 zz = 0. +1.*a;
2401 r = sqrt((xx-x)*(xx-x)+(yy-y)*(yy-y)+(zz-z)*(zz-z));
2402 return r;
2403}
a
void a()
Definition: check_aligned.C:59
s
s
Definition: check_shower.C:55
r
void r(int rid=2)
Definition: test.C:201

◆ CalcIP() [2/2]

double EdbShowerAlg_GS::CalcIP ( EdbSegP *  s,
EdbVertex *  v 
)

calculate IP between a given segment and a given vertex.
return the IP value.
this is used for IP cut.

if vertex is not given, use the selected vertex.
if(v==NULL) v=gEDA->GetSelectedVertex();

2407 {
2411
2414
2415 if (v==NULL) return -1.;
2416 return CalcIP(s, v->X(), v->Y(), v->Z());
2417}
EdbShowerAlg_GS::CalcIP
Double_t CalcIP(EdbSegP *s, double x, double y, double z)
Definition: EdbShowerAlg.cxx:2387
NULL
#define NULL
Definition: nidaqmx.h:84

◆ CheckCleanPairs()

TObjArray * EdbShowerAlg_GS::CheckCleanPairs ( EdbSegP *  InBT,
TObjArray *  RecoShowerArrayFromFindPairs 
)

Check Clean Pairs Function.
It looks for basetrack pairs that come in a row directly behind each other.
Then it will take the one most upstream in this line.
More precise to explain:: WHAT DOES THIS FUNCTION DO ???

2185{
2190
2191 Log(3,"EdbShowerAlg_GS::CheckCleanPairs","Starting CheckCleanPairs() now...");
2192
2193 if (NULL==RecoShowerArray) return NULL;
2194
2195 TObjArray* NewRecoShowerArray= new TObjArray();
2196 NewRecoShowerArray->Clear();
2197 int NewRecoShowerArrayN=0;
2198
2199 EdbTrackP* TrackPair1=NULL;
2200 EdbTrackP* TrackPair2=NULL;
2201
2202 int ntrack=RecoShowerArray->GetEntries();
2203 // cout << "EdbShowerAlg_GS::CheckCleanPairs ntrack = " << ntrack << endl;
2204
2205 for (Int_t pat_one_cnt=0; pat_one_cnt<ntrack; ++pat_one_cnt) {
2206 if (gEDBDEBUGLEVEL>2) cout << "CheckCleanPairs Doing pat_one_cnt = " << pat_one_cnt << endl;
2207 // (if -1) then the last one is not checked
2208 TrackPair1=(EdbTrackP*)RecoShowerArray->At(pat_one_cnt);
2209 bool taketrack1=true;
2210 for (Int_t pat_two_cnt=pat_one_cnt; pat_two_cnt<ntrack; ++pat_two_cnt) {
2211 if (gEDBDEBUGLEVEL>3) cout << "CheckCleanPairs Doing pat_two_cnt = " << pat_two_cnt << endl;
2212 //if only one track at all take it anyway:
2213 if (ntrack==1) continue;
2214 TrackPair2=(EdbTrackP*)RecoShowerArray->At(pat_two_cnt);
2215
2216 // Check if track1 has a track before (smaller Z value)
2217 // and if so, if dtheta is smaller than 0.1:
2218 // If both is so, then track 1 is NOT taken for final array:
2219 // cout << TrackPair1->Z() << " " << TrackPair2->Z() << endl;
2220 EdbSegP* s1=(EdbSegP* )TrackPair1->GetSegment(0);
2221 EdbSegP* s2=(EdbSegP* )TrackPair2->GetSegment(0);
2222 if (TrackPair1->Z()>TrackPair2->Z() && DeltaThetaSingleAngles((EdbSegP*)s1,(EdbSegP*)s2)<0.1) taketrack1=false;
2223 if (!taketrack1) break;
2224 }
2225
2226 if (!taketrack1) continue;
2227
2228 // Add TrackPair1
2229 NewRecoShowerArray->Add(TrackPair1);
2230 ++NewRecoShowerArrayN;
2231
2232 // Print
2233 if (TrackPair1->N()<2) cout << "CheckCleanPairs This Track has only ONE entry" << endl;
2234
2235 }
2236 if (gEDBDEBUGLEVEL>2) cout << "CheckCleanPairs From " << ntrack << " originally, there are now after Zposition and overlap cuts: " << NewRecoShowerArrayN << " left." << endl;
2237
2238 Log(3,"EdbShowerAlg_GS::CheckCleanPairs","Starting CheckCleanPairs() now...done.");
2239 return NewRecoShowerArray;
2240}
RecoShowerArray
TObjArray * RecoShowerArray
Definition: Shower_E_FromShowerRoot.C:12
EdbSegP
Definition: EdbSegP.h:21
EdbSegP::Z
Float_t Z() const
Definition: EdbSegP.h:153
EdbShowerAlg::DeltaThetaSingleAngles
Double_t DeltaThetaSingleAngles(EdbSegP *s1, EdbSegP *s2)
Definition: EdbShowerAlg.cxx:625
EdbTrackP
Definition: EdbPattern.h:113
EdbTrackP::N
Int_t N() const
Definition: EdbPattern.h:177
EdbTrackP::GetSegment
EdbSegP * GetSegment(int i) const
Definition: EdbPattern.h:195
s1
EdbSegP * s1
Definition: tlg2couples.C:29
s2
EdbSegP * s2
Definition: tlg2couples.C:30

◆ CheckInput()

Bool_t EdbShowerAlg_GS::CheckInput ( )

Before execution of the main reco routine, we perform a
check, if either an Initiator Vertex (array) was set,
or if an Initiator Basetrack (array) was set.
In this case, return kTRUE, else return kFALSE.

1246{
1251
1252 Log(2,"EdbShowerAlg_GS::CheckInput","CheckInput()");
1253
1254 Bool_t IsInput=kFALSE;
1255 Bool_t VtxArray_Or_InBTArray=kFALSE;
1256 Log(2,"EdbShowerAlg_GS::CheckInput","CheckInput() Check for eInBTArray:");
1257
1258
1259 if (eInBTArrayN==0) {
1260 Log(2,"EdbShowerAlg_GS::CheckInput","CheckInput() No eInBTArray. Check for eInVtxArray:");
1261
1262 // Check if there is an Convert_InVtxArray_To_InBTArray();
1263 if (eInVtxArrayN==0) {
1264 Log(2,"EdbShowerAlg_GS::CheckInput","CheckInput() No eInVtxArray.");
1265 Log(2,"EdbShowerAlg_GS::CheckInput","CheckInput() Set ALL BTs from the volume as InBTs:");
1266 Log(2,"EdbShowerAlg_GS::CheckInput","CheckInput() Set BTs from the volume as InBTs:");
1267 if (eFindPairsPreselected==kTRUE) {
1268 Log(2,"EdbShowerAlg_GS::CheckInput","CheckInput() Take Preselected Pairs found from");
1269 Log(2,"EdbShowerAlg_GS::CheckInput","CheckInput() the function FindPairsPreselected():");
1270 }
1271 VtxArray_Or_InBTArray=kFALSE;
1272 // Convert_EdbPVRec_To_InBTArray();
1273 // Taking all BTs takes very long time, we rather skip it.
1274 //______________________________________________________________________________
1275 // New we try finding already preselcted InBTs.... might speed up the style...
1276 if (eFindPairsPreselected==kTRUE) {
1277 SetInBTArray(this->FindPairsPreselected(eAli));
1278 }
1279 else {
1280 Convert_EdbPVRec_To_InBTArray();
1281 }
1282 eInBTArrayN=eInBTArray->GetEntries();
1283 Log(2,"EdbShowerAlg_GS::CheckInput","CheckInput() eInBTArray->GetEntries() = %d .",eInBTArray->GetEntries());
1284 //______________________________________________________________________________
1285 Log(2,"EdbShowerAlg_GS::CheckInput","CheckInput() FindPairsPreselected() done.");
1286 //eRecoMode=1;
1287// eRecoMode=2;
1288 Log(2,"EdbShowerAlg_GS::CheckInput","CheckInput() Use default eRecoMode=%d",eRecoMode);
1289 Log(2,"EdbShowerAlg_GS::CheckInput","CheckInput() Do pair search now, but take care, this is gonna take _______VERY LONG_______ time.");
1290 }
1291 else {
1292 Log(2,"EdbShowerAlg_GS::CheckInput","CheckInput() eInVtxArray there. Converting to InBTArray() now...");
1293 Convert_InVtxArray_To_InBTArray();
1294 VtxArray_Or_InBTArray=kTRUE;
1295 eRecoMode=0;
1296 Log(2,"EdbShowerAlg_GS::CheckInput","CheckInput() eRecoMode=0");
1297 Log(2,"EdbShowerAlg_GS::CheckInput","CheckInput() eInVtxArray there. Converting to InBTArray() now....done.");
1298 }
1299 }
1300
1301 cout << "EdbShowerAlg_GS::CheckInput: Finally take eRecoMode = " << eRecoMode << endl;
1302
1303 // Reload ANN weight file:
1304 // If vertices are given we will take the trained neural net for basetrack pairs
1305 // with the known vertex.
1306 // If no vertices are given then we have to assume that we dont know
1307 // what the vertex position yet is and then we look for basetrack pairings
1308 // by themselves. Therefore we need the ANN weight file trained on no vertex info.
1309
1310 ReloadANNs(eRecoMode);
1311
1312
1313 // For the MC Training, we use only highest P basetracks, if they come from signal
1314 TObjArray* arr = SelectHighestPInMCArray(eInBTArray);
1315// SetInBTArray(arr);
1316 // For the MC Training, we use only highest P basetracks, if they come from signal
1317
1318
1319 Log(2,"EdbShowerAlg_GS::CheckInput","CheckInput()...done.");
1320 return IsInput;
1321}
EdbShowerAlg_GS::SelectHighestPInMCArray
TObjArray * SelectHighestPInMCArray(TObjArray *BTArray)
Definition: EdbShowerAlg.cxx:1174
EdbShowerAlg_GS::Convert_InVtxArray_To_InBTArray
void Convert_InVtxArray_To_InBTArray()
Definition: EdbShowerAlg.cxx:1067
EdbShowerAlg_GS::eRecoMode
Int_t eRecoMode
Definition: EdbShowerAlg.h:204
EdbShowerAlg_GS::eFindPairsPreselected
Bool_t eFindPairsPreselected
Definition: EdbShowerAlg.h:201
EdbShowerAlg_GS::ReloadANNs
void ReloadANNs(Int_t RecoMode)
Definition: EdbShowerAlg.cxx:993
EdbShowerAlg_GS::FindPairsPreselected
Bool_t FindPairsPreselected(EdbSegP *InitiatorBT, EdbPVRec *eAli_Sub)
Definition: EdbShowerAlg.cxx:1674
EdbShowerAlg::SetInBTArray
void SetInBTArray(TObjArray *InBTArray)
Definition: EdbShowerAlg.h:99
EdbShowerAlg::eAli
EdbPVRec * eAli
Definition: EdbShowerAlg.h:48
EdbShowerAlg::eInBTArrayN
Int_t eInBTArrayN
Definition: EdbShowerAlg.h:52
EdbShowerAlg::Convert_EdbPVRec_To_InBTArray
void Convert_EdbPVRec_To_InBTArray()
Definition: EdbShowerAlg.cxx:1116
EdbShowerAlg::eInBTArray
TObjArray * eInBTArray
Definition: EdbShowerAlg.h:51

◆ CheckPairDuplications()

Bool_t EdbShowerAlg_GS::CheckPairDuplications ( Int_t  SegPID,
Int_t  SegID,
Int_t  Seg2PID,
Int_t  Seg2ID,
TArrayI *  SegmentPIDArray,
TArrayI *  SegmentIDArray,
TArrayI *  Segment2PIDArray,
TArrayI *  Segment2IDArray,
Int_t  RecoShowerArrayN 
)
1159{
1160 Log(3,"EdbShowerAlg_GS::CheckPairDuplications","CheckPairDuplications()");
1161 if (gEDBDEBUGLEVEL>2) cout << "EdbShowerAlg_GS::CheckPairDuplications for " << SegPID << "," << SegID <<","<< Seg2PID <<"," << Seg2ID << " compare with " << RecoShowerArrayN << " pairs:" << endl;
1162 for (Int_t i=0; i<RecoShowerArrayN; i++) {
1163 // PID and ID of Seg and Seg2 to be exchanged for duplications
1164 if ( SegPID==eSegment2PIDArray->At(i) && Seg2PID==eSegmentPIDArray->At(i) && SegID==eSegment2IDArray->At(i) && Seg2ID==eSegmentIDArray->At(i)) {
1165 if (gEDBDEBUGLEVEL>2) cout << "EdbShowerAlg_GS::CheckPairDuplications Found duplication for ... return true"<<endl;
1166 return kTRUE;
1167 }
1168 }
1169 return kFALSE;
1170}
EdbShowerAlg_GS::eSegment2PIDArray
TArrayI * eSegment2PIDArray
Definition: EdbShowerAlg.h:236
EdbShowerAlg_GS::eSegment2IDArray
TArrayI * eSegment2IDArray
Definition: EdbShowerAlg.h:237
EdbShowerAlg_GS::eSegmentPIDArray
TArrayI * eSegmentPIDArray
Definition: EdbShowerAlg.h:234
EdbShowerAlg_GS::eSegmentIDArray
TArrayI * eSegmentIDArray
Definition: EdbShowerAlg.h:235

◆ ClassDef()

EdbShowerAlg_GS::ClassDef ( EdbShowerAlg_GS  ,
1   
)

◆ Convert_InVtxArray_To_InBTArray()

void EdbShowerAlg_GS::Convert_InVtxArray_To_InBTArray ( )

This function takes vertices from the InVtxArray object,
and fills the InBTArray with those: X/Y/Z position as
from the Vertices, TX(),TY() equal to 0
(this doesnt matter when having InVtxArray with vertices filled
for this algo, since it does rely only on the point
position here).
This function is useful, if one wants to start from
already calculated vertices.
If vtx-array is not filled, on the other hand, it is recommended
to use the full EdbPVRec object as possible starting points,
look in EdbShowerAlg::Convert_EdbPVRec_To_InBTArray().

1068{
1081
1082 Log(2,"EdbShowerAlg_GS::EdbShowerAlg_GS","Convert_InVtxArray_To_InBTArray()");
1083 if (eInBTArrayN!=0) {
1084 eInBTArray->Clear();
1085 cout << " eInBTArray->Clear();" << endl;
1086 }
1087
1088 if (eInVtxArray==NULL || eInVtxArrayN==0 || eInVtxArraySet==kFALSE) cout << "NO eInVtxArray " << endl;
1089 EdbVertex* vtx;
1090
1091 Log(2,"EdbShowerAlg_GS::EdbShowerAlg_GS","Convert_InVtxArray_To_InBTArray()...start loop (N=%d):",eInVtxArrayN);
1092 for (Int_t i=0; i<eInVtxArrayN; i++) {
1093 vtx= (EdbVertex*)eInVtxArray->At(i);
1094 EdbSegP* seg = new EdbSegP(i,vtx->X(),vtx->Y(),0,0);
1095 seg->SetZ(vtx->Z());
1096 // vtx can have by default initialization MCEvt==0,
1097 // this we dont want, in case: we set MCEvt of vtx from 0 to -999
1098 if (vtx->MCEvt()==0) vtx->SetMC(-999);
1099 seg->SetMC(vtx->MCEvt(),vtx->MCEvt());
1100 seg->SetFlag(0);
1101 eInBTArray->Add(seg);
1102
1103 //cout <<"DEBUG vtx->MCEvt() " << vtx->MCEvt() << endl;
1104 // seg->PrintNice();
1105 }
1106
1107 eInBTArrayN=eInBTArray->GetEntries();
1108 if (gEDBDEBUGLEVEL>2) cout << "EdbShowerAlg_GS::Convert_InVtxArray_To_InBTArray Converted " << eInBTArrayN << "InVtx to InBT." << endl;
1109 Log(2,"EdbShowerAlg_GS::EdbShowerAlg_GS","Convert_InVtxArray_To_InBTArray()...done.");
1110 return;
1111}
EdbSegP::SetZ
void SetZ(float z)
Definition: EdbSegP.h:125
EdbSegP::SetMC
void SetMC(int mEvt, int mTrack)
Definition: EdbSegP.h:141
EdbSegP::SetFlag
void SetFlag(int flag)
Definition: EdbSegP.h:130
EdbVertex
Definition: EdbVertex.h:69
EdbVertex::SetMC
void SetMC(int mEvt=0)
Definition: EdbVertex.h:159

◆ CreateANNPair()

void EdbShowerAlg_GS::CreateANNPair ( )

Create the neural network, which is used as last layer in the
GS algo reconstruction.
Load its weights from the $FEDRASYS path.
Initial testings of the method have been done using the ROOT
TMVA framework, espacially the multivariate method "Boosted Decision Trees"
and the two given Neural Network Implementations and linear cuts.
ROC curves can be found in the manual. The two NN implementation dont
differ too much, and we can take -with good conscience- the standard
ROOT TMultiLayerPerceptron one.

903{
913
914 Log(2,"EdbShowerAlg_GS::EdbShowerAlg_GS","CreateANNPair()");
915
916 // Just to avoid the root errormessages of creating a tree non memory resident.
917 TFile* localANNPairTreeFile = new TFile("localANNPairTreeFile.root","RECREATE");
918
919 // Create the tree which will store the Inputvariables for the GS ANN discrimination:
920 eANNPairTree= new TTree("eANNPairTree","eANNPairTree");
921 cout << "EdbShowerAlg_GS::CreateANNPair() eANNPairTree Memory creation done." << endl;
922 eANNPairTree->Branch("eValueGSNN_varInput",&eValueGSNN_varInput,"eValueGSNN_varInput/F");
923 eANNPairTree->Branch("eValueGSNN_varOutput",&eValueGSNN_varOutput,"eValueGSNN_varOutput/F");
924 eANNPairTree->Branch("eValueGSNN_var00",&eValueGSNN_var00,"eValueGSNN_var00/F");
925 eANNPairTree->Branch("eValueGSNN_var01",&eValueGSNN_var01,"eValueGSNN_var01/F");
926 eANNPairTree->Branch("eValueGSNN_var02",&eValueGSNN_var02,"eValueGSNN_var02/F");
927 eANNPairTree->Branch("eValueGSNN_var03",&eValueGSNN_var03,"eValueGSNN_var03/F");
928 eANNPairTree->Branch("eValueGSNN_var04",&eValueGSNN_var04,"eValueGSNN_var04/F");
929 eANNPairTree->Branch("eValueGSNN_var05",&eValueGSNN_var05,"eValueGSNN_var05/F");
930 eANNPairTree->Branch("eValueGSNN_var06",&eValueGSNN_var06,"eValueGSNN_var06/F");
931 cout << "EdbShowerAlg_GS::CreateANNPair() eANNPairTree SetBranchAddress done." << endl << endl;
932
933 // :@eValueGSNN_var06 is not used, because Flag() is a MC information only.
934 TString layout="";
935 TString weights="";
936
937 // Attention: the three cases have a different number of input neurons:
938 // Case A and Case C have 6 input neurons
939 // Case B has only 5 input neurons (because dZ is a constant and we cannot use ANN
940 // with constant input neurons).
941 //
942 // CASE A : Vertex Info out of "RealVertex"
943 // CASE B : No Vertex Info, but "FakeVertex" made out of one BT, dZ==650
944 // CASE C : No Vertex Info, but "InterimVertex" made out of two BTs, dZ variable
945
946 // CASE A :
947 layout="@eValueGSNN_var00,@eValueGSNN_var01,@eValueGSNN_var02,@eValueGSNN_var03,@eValueGSNN_var04,@eValueGSNN_var05:7:6:eValueGSNN_varInput";
948 // Create the Multilayerperceptron
949 eANNPairCaseA = new TMultiLayerPerceptron(layout,eANNPairTree,"","");
950 cout << "EdbShowerAlg_GS::CreateANNPair() eANNPairCaseA TMultiLayerPerceptron creation done." << endl;
951 cout << "EdbShowerAlg_GS::CreateANNPair() Layout = " << endl;
952 cout << layout.Data() << endl;
953 weights= TString(gSystem->ExpandPathName("$FEDRA_ROOT"))+TString("/src/libShower/weights/Reco/ShowerAlg_GS/weights_CASE_A.txt");
954 eANNPairCaseA->LoadWeights(weights);
955 cout << "EdbShowerAlg_GS::CreateANNPair() eANNPair: weights = " << endl;
956 cout << weights.Data() << endl << endl;
957
958 // CASE B :
959 layout="@eValueGSNN_var01,@eValueGSNN_var02,@eValueGSNN_var04,@eValueGSNN_var05:7:6:eValueGSNN_varInput";
960 // Create the Multilayerperceptron
961 eANNPairCaseB = new TMultiLayerPerceptron(layout,eANNPairTree,"","");
962 cout << "EdbShowerAlg_GS::CreateANNPair() eANNPairCaseB TMultiLayerPerceptron creation done." << endl;
963 cout << "EdbShowerAlg_GS::CreateANNPair() Layout = " << endl;
964 cout << layout.Data() << endl;
965 weights= TString(gSystem->ExpandPathName("$FEDRA_ROOT"))+TString("/src/libShower/weights/Reco/ShowerAlg_GS/weights_CASE_B.txt");
966 eANNPairCaseB->LoadWeights(weights);
967 cout << "EdbShowerAlg_GS::CreateANNPair() eANNPair: weights = " << endl;
968 cout << weights.Data() << endl << endl;
969
970 // CASE C :
971 layout="@eValueGSNN_var00,@eValueGSNN_var01,@eValueGSNN_var02,@eValueGSNN_var03,@eValueGSNN_var04,@eValueGSNN_var05:7:6:eValueGSNN_varInput";
972 // Create the Multilayerperceptron
973 eANNPairCaseC = new TMultiLayerPerceptron(layout,eANNPairTree,"","");
974 cout << "EdbShowerAlg_GS::CreateANNPair() eANNPairCaseC TMultiLayerPerceptron creation done." << endl;
975 cout << "EdbShowerAlg_GS::CreateANNPair() Layout = " << endl;
976 cout << layout.Data() << endl;
977 weights= TString(gSystem->ExpandPathName("$FEDRA_ROOT"))+TString("/src/libShower/weights/Reco/ShowerAlg_GS/weights_CASE_C.txt");
978 eANNPairCaseC->LoadWeights(weights);
979 cout << "EdbShowerAlg_GS::CreateANNPair() eANNPair: weights = " << endl;
980 cout << weights.Data() << endl << endl;
981
982 // By default we take the weights for the Case C (explained in the manual)
983 // since we dont know a priori, if vertices are given.
984 eANNPair=eANNPairCaseC;
985
986 Log(2,"EdbShowerAlg_GS::EdbShowerAlg_GS","CreateANNPair()...done.");
987 return;
988}
EdbShowerAlg_GS::eValueGSNN_var02
Float_t eValueGSNN_var02
Definition: EdbShowerAlg.h:225
EdbShowerAlg_GS::eANNPairCaseC
TMultiLayerPerceptron * eANNPairCaseC
Definition: EdbShowerAlg.h:219
EdbShowerAlg_GS::eANNPairCaseB
TMultiLayerPerceptron * eANNPairCaseB
Definition: EdbShowerAlg.h:218
EdbShowerAlg_GS::eANNPair
TMultiLayerPerceptron * eANNPair
Definition: EdbShowerAlg.h:216
EdbShowerAlg_GS::eValueGSNN_var05
Float_t eValueGSNN_var05
Definition: EdbShowerAlg.h:228
EdbShowerAlg_GS::eValueGSNN_var03
Float_t eValueGSNN_var03
Definition: EdbShowerAlg.h:226
EdbShowerAlg_GS::eANNPairTree
TTree * eANNPairTree
Definition: EdbShowerAlg.h:220
EdbShowerAlg_GS::eValueGSNN_var04
Float_t eValueGSNN_var04
Definition: EdbShowerAlg.h:227
EdbShowerAlg_GS::eValueGSNN_var01
Float_t eValueGSNN_var01
Definition: EdbShowerAlg.h:224
EdbShowerAlg_GS::eValueGSNN_var00
Float_t eValueGSNN_var00
Definition: EdbShowerAlg.h:223
EdbShowerAlg_GS::eValueGSNN_varOutput
Float_t eValueGSNN_varOutput
Definition: EdbShowerAlg.h:231
EdbShowerAlg_GS::eANNPairCaseA
TMultiLayerPerceptron * eANNPairCaseA
Definition: EdbShowerAlg.h:217
EdbShowerAlg_GS::eValueGSNN_var06
Float_t eValueGSNN_var06
Definition: EdbShowerAlg.h:229
EdbShowerAlg_GS::eValueGSNN_varInput
Float_t eValueGSNN_varInput
Definition: EdbShowerAlg.h:230

◆ CreateANNPlots()

void EdbShowerAlg_GS::CreateANNPlots ( )

Just a plot creating function for the GS_ANN Input variables.

2343{
2345
2346 Log(3,"EdbShowerAlg_GS::CreateANNPlots","CreateANNPlots...");
2347
2348 gROOT->SetStyle("Bold");
2349 // Draw Histograms for the last cut on Gamma Search
2350 TCanvas* canv_input = new TCanvas();
2351 canv_input->Divide(3,2);
2352 canv_input->cd(1);
2353 eANNPairTree->Draw("eValueGSNN_var00");
2354 canv_input->cd(2);
2355 eANNPairTree->Draw("eValueGSNN_var01");
2356 canv_input->cd(3);
2357 eANNPairTree->Draw("eValueGSNN_var02");
2358 canv_input->cd(4);
2359 eANNPairTree->Draw("eValueGSNN_var03");
2360 canv_input->cd(5);
2361 eANNPairTree->Draw("eValueGSNN_var04");
2362 canv_input->cd(6);
2363 eANNPairTree->Draw("eValueGSNN_var05");
2364
2365 TCanvas* canv_InOutput = new TCanvas();
2366 canv_InOutput->Divide(3,1);
2367 canv_InOutput->cd(1);
2368 eANNPairTree->Draw("eValueGSNN_varInput");
2369 canv_InOutput->cd(2);
2370 eANNPairTree->Draw("eValueGSNN_varOutput");
2371 canv_InOutput->cd(3);
2372 eANNPairTree->Draw("eValueGSNN_varInput:eValueGSNN_varOutput","","colz");
2373
2374 Log(3,"EdbShowerAlg_GS::CreateANNPlots","CreateANNPlots...done.");
2375 return;
2376}
TCanvas
new TCanvas()

◆ Execute()

void EdbShowerAlg_GS::Execute ( )
virtual

Main shower reco function.
Loops over all filled Initiator basetracks.
Fills the TObjArray with reconstructed showers.
GS Reco algo is in function FindPairs().

Debug Purpose:

Reimplemented from EdbShowerAlg.

1379{
1384
1385 Log(2,"EdbShowerAlg_GS::Execute","Execute()");
1386 Log(2,"EdbShowerAlg_GS::Execute","Execute() DOING MAIN SHOWER RECONSTRUCTION HERE");
1387 Log(2,"EdbShowerAlg_GS::Execute","Execute() Check for eInBTArray:");
1388
1389 // Check Input which BT filling should be done.
1390 Bool_t VtxArray_Or_InBTArray = kFALSE;
1391 VtxArray_Or_InBTArray = CheckInput();
1392
1393
1394 //--- Needed interim objects:
1395 EdbSegP* InBT=NULL;
1396// EdbShowerP* RecoShower; //libShowRec
1397 EdbTrackP* RecoShower; //libShower
1398
1399 Int_t PermilleCount=eInBTArrayN/1000;
1400 if (PermilleCount==0) PermilleCount=1; // floating point exception failsafe.
1401
1402 Int_t STEP=-1;
1403 if (eFirstPlate_eAliPID-eLastPlate_eAliPID<0) STEP=1;
1404 if (gEDBDEBUGLEVEL>2) cout << "EdbShowerAlg_GS::Execute--- STEP for patternloop direction = " << STEP << endl;
1405 if (gEDBDEBUGLEVEL>2) cout << "EdbShowerAlg_GS::Execute--- eRecoShowerArrayN = " << eRecoShowerArrayN << endl;
1406
1407 //--- Loop over InBTs:
1408 if (gEDBDEBUGLEVEL>=2) {
1409 cout << "EdbShowerAlg_GS::Execute Loop over InBTs: " << endl;
1410 cout << "EdbShowerAlg_GS::Execute (N=" << eInBTArrayN << "):" << endl;
1411 cout << "EdbShowerAlg_GS::Execute (each dot means a progress of 1 permille of N):" << endl;
1412 cout << "EdbShowerAlg_GS::Execute (each K means a progress of 1 K InBTs):" << endl;
1413 }
1414
1415 //-----------------------------------
1416 // Since eInBTArray is filled in ascending ordering by zpositon
1417 // We use the descending loop to begin with BT with lowest z first.
1418 // InBT doest have to be necessary a real BaseTrack, it can also be a vertex (with its positions and angle zero) !!!
1419 for (Int_t i=eInBTArrayN-1; i>=0; --i) {
1420
1421 if (eDebug) if (i<eInBTArrayN-1) continue;
1422
1423 //-----------------------------------
1424 // 0) Set counters...
1425 //-----------------------------------
1426 Log(3,"EdbShowerAlg_GS::Execute","//-----------------------------------");
1427 Log(3,"EdbShowerAlg_GS::Execute","// 0) Set counters...:");
1428 Log(3,"EdbShowerAlg_GS::Execute","//-----------------------------------");
1429
1430 // CounterOutPut
1431// if (gEDBDEBUGLEVEL==2) if ((i%100)==0) cout << eInBTArrayN <<" InBT in total, still to do:"<<Form("%4d",i)<< "\r\r\r\r"<<flush;
1432 if (gEDBDEBUGLEVEL==2) if ((i%PermilleCount)==0) cout << "." << flush;
1433 if (gEDBDEBUGLEVEL==2) if ((i%1000)==0) cout << "K" << flush;
1434
1435 // Get InitiatorBT from eInBTArray InBT
1436 InBT=(EdbSegP*)eInBTArray->At(i);
1437 if (gEDBDEBUGLEVEL>2) {
1438 Log(3,"EdbShowerAlg_GS::Execute","//---- Print Initiator BT(Vtx) ------------------");
1439 InBT->PrintNice();
1440 }
1441
1442 //-----------------------------------
1443 // 1) Make local_gAli with cut parameters:
1444 //-----------------------------------
1445 Log(3,"EdbShowerAlg_GS::Execute","//-----------------------------------");
1446 Log(3,"EdbShowerAlg_GS::Execute","// 1) Make local_gAli with cut parameters:");
1447 Log(3,"EdbShowerAlg_GS::Execute","//-----------------------------------");
1448 // Transform (make size smaller, extract only events having same MC) the eAli object:
1449 Transform_eAli(InBT,999999);
1450 if (gEDBDEBUGLEVEL>2) eAli_Sub->Print();
1451
1452
1453 //-----------------------------------
1454 // 2) FindPairs
1455 //-----------------------------------
1456 Log(3,"EdbShowerAlg_GS::Execute","//-----------------------------------");
1457 Log(3,"EdbShowerAlg_GS::Execute","// 2) FindPairs");
1458 Log(3,"EdbShowerAlg_GS::Execute","//-----------------------------------");
1459 TObjArray* Pairs = FindPairs(InBT,eAli_Sub);
1460 if (gEDBDEBUGLEVEL>2) cout << "TObjArray* Pairs = FindPairs(InBT,eAli_Sub); done. Entries= " << Pairs->GetEntries() << endl;
1461
1462
1463
1464 //-----------------------------------
1465 // 3) Clear Found Pairs
1466 //-----------------------------------
1467 Log(3,"EdbShowerAlg_GS::Execute","//-----------------------------------");
1468 Log(3,"EdbShowerAlg_GS::Execute","// 3) Clear Found Pairs");
1469 Log(3,"EdbShowerAlg_GS::Execute","//-----------------------------------");
1470 TObjArray* CleanPairs;
1471 if (eSetCleanPairs==kTRUE) {
1472 CleanPairs = CheckCleanPairs( InBT, Pairs );
1473 // After we found pairs (after clear) we add the internal RecoShowerArray:
1474 AddRecoShowerArray(CleanPairs);
1475 FillANNTree(CleanPairs,InBT);
1476 }
1477 else {
1478 // After we found pairs (without clear) we add the internal RecoShowerArray:
1479 AddRecoShowerArray(Pairs);
1480 FillANNTree(Pairs,InBT);
1481 }
1482
1483 //-----------------------------------
1484 // 4) Finally Apply GS_ANN cut on Pairs
1485 // ---> Is already done in the FindPairs() routine!!
1486 //-----------------------------------
1487
1488
1489 } // Loop over InBT
1490 //-----------------------------------
1491
1492 cout << "EdbShowerAlg_GS::Execute Loop over InBT finished."<< endl;
1493 cout << "EdbShowerAlg_GS::Execute eRecoShowerArray= " << eRecoShowerArray << endl;
1494 cout << "EdbShowerAlg_GS::Execute eRecoShowerArrayN= " << eRecoShowerArrayN << endl;
1495
1496
1497 if (gEDBDEBUGLEVEL>2) {
1498 cout << "EdbShowerAlg_GS::Execute Print all BT pairs now:" << endl;
1499 for (Int_t i=0; i <eRecoShowerArray->GetEntries(); i++) {
1500// EdbShowerP* sh = (EdbShowerP* )eRecoShowerArray->At(i); // libShowRec
1501 EdbTrackP* sh = (EdbTrackP* )eRecoShowerArray->At(i); // libShower
1502 sh->PrintNice();
1503 }
1504 }
1505
1507 if (eRecoShowerArrayN>0) {
1508 EdbTrackP* sh = (EdbTrackP* )eRecoShowerArray->At(0); // libShower
1509 cout << "EdbShowerAlg_GS::Execute Shower 0: PrintNice: " << endl;
1510 sh->PrintNice();
1511 sh = (EdbTrackP* )eRecoShowerArray->At(eRecoShowerArray->GetEntries()-1); // libShower
1512 cout << "EdbShowerAlg_GS::Execute Shower last("<< eRecoShowerArrayN <<"): PrintNice: " << endl;
1513 sh->PrintNice();
1514 }
1515
1516 eANNPair->Write();
1517 eANNPairTree->Write();
1518
1519 Log(2,"EdbShowerAlg_GS::Execute","Execute()...done.");
1520 return;
1521}
EdbPatternsVolume::Print
void Print() const
Definition: EdbPattern.cxx:1693
EdbSegP::PrintNice
void PrintNice() const
Definition: EdbSegP.cxx:392
EdbShowerAlg_GS::eSetCleanPairs
Bool_t eSetCleanPairs
Definition: EdbShowerAlg.h:210
EdbShowerAlg_GS::CheckCleanPairs
TObjArray * CheckCleanPairs(EdbSegP *InBT, TObjArray *RecoShowerArrayFromFindPairs)
Definition: EdbShowerAlg.cxx:2184
EdbShowerAlg_GS::CheckInput
Bool_t CheckInput()
Definition: EdbShowerAlg.cxx:1245
EdbShowerAlg_GS::FindPairs
TObjArray * FindPairs(EdbSegP *InBT, EdbPVRec *eAli_Sub)
Definition: EdbShowerAlg.cxx:1751
EdbShowerAlg_GS::FillANNTree
void FillANNTree(TObjArray *RecoShowerArray, EdbSegP *Inbt)
Definition: EdbShowerAlg.cxx:2245
EdbShowerAlg::eRecoShowerArray
TObjArray * eRecoShowerArray
Definition: EdbShowerAlg.h:59
EdbShowerAlg::AddRecoShowerArray
void AddRecoShowerArray(TObjArray *RecoShowerArray)
Definition: EdbShowerAlg.cxx:108
EdbShowerAlg::eLastPlate_eAliPID
Int_t eLastPlate_eAliPID
Definition: EdbShowerAlg.h:55
EdbShowerAlg::Transform_eAli
void Transform_eAli(EdbSegP *InitiatorBT, Float_t ExtractSize)
Definition: EdbShowerAlg.cxx:137
EdbShowerAlg::eFirstPlate_eAliPID
Int_t eFirstPlate_eAliPID
Definition: EdbShowerAlg.h:54
EdbShowerAlg::eDebug
Bool_t eDebug
Definition: EdbShowerAlg.h:80
EdbShowerAlg::eRecoShowerArrayN
Int_t eRecoShowerArrayN
Definition: EdbShowerAlg.h:60
EdbShowerAlg::eAli_Sub
EdbPVRec * eAli_Sub
Definition: EdbShowerAlg.h:63
EdbTrackP::PrintNice
void PrintNice()
Definition: EdbPattern.cxx:1169

◆ FillANNTree()

void EdbShowerAlg_GS::FillANNTree ( TObjArray *  RecoShowerArray,
EdbSegP *  Inbt 
)

only here we distinguish, since the function DeltaR_WithPropagation is sensitive to the Z-ordering which basetrack comes first.

2246{
2247 // This Function will be used to fill the ANN tree,
2248 // when the routine is called _after_ FindPairs.
2249 // Value calulation is the same as in the FindPairs function.
2250
2251 Log(3,"EdbShowerAlg_GS::FillANNTree","FillANNTree...");
2252 Log(3,"EdbShowerAlg_GS::FillANNTree","RecoShowerArray->GetEntries() = %d", RecoShowerArray->GetEntries());
2253
2254 if (RecoShowerArray->GetEntries()==0) return;
2255
2256 for (Int_t i=0; i< RecoShowerArray->GetEntries(); ++i) {
2257
2258 EdbTrackP* shower = (EdbTrackP*)RecoShowerArray->At(i);
2259 EdbSegP* Segment = (EdbSegP*)shower->GetSegment(0);
2260 EdbSegP* Segment2 = (EdbSegP*)shower->GetSegment(1);
2261
2262 EdbVertex* vtx=new EdbVertex();
2263 if (eRecoMode==1 || eRecoMode==2 ) {
2264 vtx->SetXYZ(InBT->X()-650*InBT->TX(),InBT->Y()-650*InBT->TY(),InBT->Z()-650);
2265 }
2266 else {
2267 vtx->SetXYZ(InBT->X(),InBT->Y(),InBT->Z());
2268 }
2269
2270 TObjArray *segments = new TObjArray(2);
2271 segments->Add(Segment);
2272 segments->Add(Segment2);
2273 EdbVertex* vetex = CalcVertex(segments);
2274 delete segments;
2275 Float_t IP_Seg1_To_Vtx = CalcIP(Segment, vtx);
2276 Float_t IP_Seg2_To_Vtx = CalcIP(Segment2,vtx);
2277 Float_t IP_InBT_To_Vtx = CalcIP(InBT,vtx);
2278 Float_t IP_Seg1ToVtxSeg1Seg2=CalcIP(Segment ,vetex);
2279 Float_t IP_Seg2ToVtxSeg1Seg2=CalcIP(Segment2,vetex);
2280
2281 eValueGSNN_varInput=-1;
2282 // Purity 1: Input = 1.0;
2283 // Purity 0.5: Input = 0.5;
2284 // Purity else: Input = 0.0;
2285 if (Segment2->Flag()+Segment->Flag()==0&&TMath::Abs(Segment2->Flag())==11&&Segment->MCEvt()>0) {
2286 eValueGSNN_varInput=1;
2287 }
2288 else if (Segment2->Flag()+Segment->Flag()!=0&&TMath::Abs(Segment2->Flag())==11&&Segment->MCEvt()>0) {
2289 eValueGSNN_varInput=0.5;
2290 }
2291 else {
2292 eValueGSNN_varInput=0;
2293 }
2294
2295 // Set Tree Input Variable Values:
2296 // Calculation of these values is EXACTLY the same as in the
2297 // FindPairs() function.
2298 eValueGSNN_var00=TMath::Min(IP_Seg1_To_Vtx,IP_Seg2_To_Vtx);
2299 if (eRecoMode==1) eValueGSNN_var00=IP_InBT_To_Vtx;
2300 if (eRecoMode==2) {
2301 eValueGSNN_var00=TMath::Min(IP_Seg1ToVtxSeg1Seg2,IP_Seg2ToVtxSeg1Seg2);
2302 }
2303 eValueGSNN_var01=GetMinimumDist(Segment,Segment2);
2304 eValueGSNN_var02=DeltaR_WithPropagation(Segment,Segment2);
2307 if (Segment2->Z()>Segment->Z()) eValueGSNN_var02=DeltaR_WithPropagation(Segment2,Segment);
2308 eValueGSNN_var03=InBT->Z()-vtx->Z();
2309 if (eRecoMode==2) {
2310 eValueGSNN_var03=TMath::Min(Segment->Z(),Segment2->Z())-vetex->Z();
2311 }
2312 if (eRecoMode==0) {
2313 eValueGSNN_var03=TMath::Min(Segment->Z(),Segment2->Z())-vtx->Z();
2314 }
2315 eValueGSNN_var04=DeltaThetaSingleAngles(Segment,Segment2);
2316 eValueGSNN_var05=TMath::Abs(Segment->PID()-Segment2->PID());
2317 eValueGSNN_var06=Segment2->Flag()+Segment->Flag();
2318
2319 if (gEDBDEBUGLEVEL>2) {
2320 cout <<"------------- FillANNTree Print eValueGSNN first, pring seg1 and seg2------------- " << endl;
2321 Segment->PrintNice();
2322 Segment2->PrintNice();
2323 cout <<"------------- FillANNTree Print eValueGSNN_var00= " << eValueGSNN_var00 << " ---- " << endl;
2324 cout <<"------------- FillANNTree Print eValueGSNN_var01= " << eValueGSNN_var01 << " ---- " << endl;
2325 cout <<"------------- FillANNTree Print eValueGSNN_var02= " << eValueGSNN_var02 << " ---- " << endl;
2326 cout <<"------------- FillANNTree Print eValueGSNN_var03= " << eValueGSNN_var03 << " ---- " << endl;
2327 cout <<"------------- FillANNTree Print eValueGSNN_var04= " << eValueGSNN_var04 << " ---- " << endl;
2328 cout <<"------------- FillANNTree Print eValueGSNN_var05= " << eValueGSNN_var05 << " ---- " << endl;
2329 }
2330
2331 eANNPairTree->Fill();
2332
2333 } // for (Int_t i=0; i< RecoShowerArray->GetEntries(); ++i) {
2334
2335 Log(3,"EdbShowerAlg_GS::FillANNTree","FillANNTree...done.");
2336 return;
2337}
EdbSegP::PID
Int_t PID() const
Definition: EdbSegP.h:148
EdbSegP::MCEvt
Int_t MCEvt() const
Definition: EdbSegP.h:145
EdbSegP::Flag
Int_t Flag() const
Definition: EdbSegP.h:149
EdbShowerAlg::DeltaR_WithPropagation
Double_t DeltaR_WithPropagation(EdbSegP *s, EdbSegP *stest)
Definition: EdbShowerAlg.cxx:577
EdbShowerAlg::GetMinimumDist
Double_t GetMinimumDist(EdbSegP *seg1, EdbSegP *seg2)
Definition: EdbShowerAlg.cxx:651
EdbShowerAlg::CalcVertex
EdbVertex * CalcVertex(TObjArray *segments)
Definition: EdbShowerAlg.cxx:353
EdbVertex::SetXYZ
void SetXYZ(float x, float y, float z)
Definition: EdbVertex.h:157

◆ Finalize()

void EdbShowerAlg_GS::Finalize ( )
virtual

do nothing yet.

Reimplemented from EdbShowerAlg.

2381{
2383}

◆ FindPairs()

TObjArray * EdbShowerAlg_GS::FindPairs ( EdbSegP *  InBT,
EdbPVRec *  eAli_Sub 
)

Function to find pairs of an EdbPVRec object.
with respect to a Initiator Basetrack (either a real vertex
or a helper vertex from the BT itself).
ATTENTION: CURRENTLY (23.07.2011) STILL IN MODIFYING PHASE!!!
ATTENTION: CURRENTLY (19.09.2011) STILL IN MODIFYING PHASE!!!
ATTENTION: CURRENTLY (01.10.2011) STILL IN MODIFYING PHASE!!!
ATTENTION: CURRENTLY (05.10.2011) STILL IN MODIFYING PHASE!!!
ATTENTION: CURRENTLY (11.10.2011) STILL IN MODIFYING PHASE!!!

if (eRecoMode==1 && distZ>10) continue; in this reco mode we loop only over segment2 , since InBT and vtx are correlated!!!

1752{
1762
1763 Log(3,"EdbShowerAlg_GS::FindPairs","Starting FindPairs(InBT,eAli_Sub) now....");
1764
1765 TObjArray* RecoShowerArray= new TObjArray();
1766 RecoShowerArray->Clear();
1767 Int_t RecoShowerArrayN=0;
1768
1769 eSegmentPIDArray->Reset();
1770 eSegmentIDArray->Reset();
1771 eSegment2PIDArray->Reset();
1772 eSegment2IDArray->Reset();
1773
1774 EdbSegP* Segment=NULL;
1775 EdbSegP* Segment2=NULL;
1776 Float_t x_av,y_av,z_av,tx_av,ty_av,distZ;
1777 EdbSegP* Segment_Sum=new EdbSegP(0,0,0,0,0,0);
1778 Float_t IP_Pair_To_InBT=0;
1779 Float_t IP_Pair_To_InBT_Seg=0;
1780 Float_t IP_Pair_To_InBT_Seg2=0;
1781 Float_t IP_Pair_To_InBT_SegSum=0;
1782 Float_t IP_Pair_To_InBT_SegSmaller=0;
1783 Float_t IP_InBT_To_Vtx=0;
1784 Float_t IP_Seg1_To_Vtx=0;
1785 Float_t IP_Seg2_To_Vtx=0;
1786
1787 Int_t NPairTriesTotal=0;
1788 Double_t paramsAC[6]; // for Case A and C
1789 Double_t paramsB[4]; // for Case B
1790
1791
1792 EdbSegP* InBT=NULL;
1793 if (NULL==InitiatorBT) {
1794 InBT= new EdbSegP();
1795 EdbPattern* pat = eAli_Sub->GetPatternZLowestHighest(1);
1796 InBT->SetX(pat->X());
1797 InBT->SetY(pat->Y());
1798 InBT->SetZ(pat->Z());
1799 InBT->SetTX(0);
1800 InBT->SetTY(0);
1801 InBT->SetMC(-999,-999);
1802 cout << "WARNING EdbShowerAlg_GS::FindPairs InBT==NULL. Created a dummy InBT:" << endl;
1803 InBT->PrintNice();
1804 }
1805 else {
1806 InBT=InitiatorBT;
1807 }
1808
1809 // Create Helper Vertex:
1810 EdbVertex* vtx=new EdbVertex();
1811 vtx->SetMC(InBT->MCEvt());
1812
1813 // Now here distinguish the Case A and Case B for the gamma reco:
1814 // Case A(0): vertex is the InBT.
1815 // Case B(1): vertex is made out of the InBT, propagated back half a plate.
1816 // Case C(2): like case B, but for the IP cutvariable, the IP to new vertex, made
1817 // out of BT 1 & 2 is used.
1818 if (eRecoMode==1 || eRecoMode==2 ) {
1819 if (gEDBDEBUGLEVEL>2) cout << "EdbShowerAlg_GS::FindPairs eRecoMode==1: Propagate InBT half a plate back for the vertex." << endl;
1820 vtx->SetXYZ(InBT->X()-650*InBT->TX(),InBT->Y()-650*InBT->TY(),InBT->Z()-650);
1821 }
1822 else {
1823 if (gEDBDEBUGLEVEL>2) cout << "EdbShowerAlg_GS::FindPairs eRecoMode==0: Take InBT position for the vertex." << endl;
1824 //cout << "DEBUG EdbShowerAlg_GS::FindPairs eRecoMode==0: Take InBT position for the vertex." << endl;
1825 vtx->SetXYZ(InBT->X(),InBT->Y(),InBT->Z());
1826 }
1827
1828 IP_InBT_To_Vtx= CalcIP(InBT,vtx);
1829 //cout << " CalcIP(InBT,vtx) " << CalcIP(InBT,vtx) << endl;
1830
1831 //-----------------------------------
1832 // 2) Loop over (whole) eAli, check BT for Cuts
1833 // eAli_Sub
1834 // Loop structure:
1835 // Loop over plates [0..NPatterns-1] for pattern one
1836 // Loop over plates [0..NPatterns-1] for pattern two
1837 // Take only plate pairings for |PID diff|<=3
1838 // Loop over all BT of pattern one
1839 // Loop over all BT of pattern two
1840
1841 // This doesnt yet distinguish the FIRSTPLATE, MIDDLEPLATE, LATPLATE, NUMBERPLATES
1842 // labelling, this will be built in later....
1843 //-----------------------------------
1844
1845 Int_t npat=eAli_Sub->Npatterns();
1846 Int_t pat_one_bt_cnt_max,pat_two_bt_cnt_max=0;
1847 EdbPattern* pat_one=0;
1848 EdbPattern* pat_two=0;
1849
1850 //-----------------------------------
1851 // --- Loop over first patterns
1852 for (Int_t pat_one_cnt=0; pat_one_cnt<npat; ++pat_one_cnt) {
1853 pat_one=(EdbPattern*)eAli_Sub->GetPattern(pat_one_cnt);
1854 pat_one_bt_cnt_max=eAli_Sub->GetPattern(pat_one_cnt)->GetN();
1855
1856 // Check if dist Z to vtx (BT) is ok:
1857 distZ=pat_one->Z()-InBT->Z();
1858 if (gEDBDEBUGLEVEL>3) {
1859 cout << "EdbShowerAlg_GS::FindPairs Check if dist Z to vtx (BT) is ok: distZ=" << distZ << endl;
1860 }
1861 if (distZ>eParaValue[3]) continue;
1862 if (distZ<-1000) continue;
1863
1866
1867 if (gEDBDEBUGLEVEL>2) cout << "EdbShowerAlg_GS::FindPairs pat_one_cnt=" << pat_one_cnt << " pat_one->Z() = " << pat_one->Z() << " pat_one_bt_cnt_max= "<< pat_one_bt_cnt_max <<endl;
1868
1869 //-----------------------------------
1870 // --- Loop over second patterns
1871 for (Int_t pat_two_cnt=0; pat_two_cnt<npat; ++pat_two_cnt) {
1872
1873 // Now apply cut conditions: GS GAMMA SEARCH Alg:
1874 if (TMath::Abs(pat_one_cnt-pat_two_cnt)>eParaValue[5]) continue;
1875
1876 pat_two=(EdbPattern*)eAli_Sub->GetPattern(pat_two_cnt);
1877 pat_two_bt_cnt_max=eAli_Sub->GetPattern(pat_two_cnt)->GetN();
1878
1879 distZ=pat_two->Z()-InBT->Z();
1880 if (distZ>eParaValue[3]) continue;
1881 if (distZ<-1000) continue;
1882
1883 if (gEDBDEBUGLEVEL>2) cout << "EdbShowerAlg_GS::FindPairs pat_two_cnt=" << pat_two_cnt << " pat_two->Z() = " << pat_two->Z() << " pat_two_bt_cnt_max= "<< pat_two_bt_cnt_max <<endl;
1884
1885 //-----------------------------------
1886 // --- Loop over basetracks of first pattern
1887 for (Int_t pat_one_bt_cnt=0; pat_one_bt_cnt<pat_one_bt_cnt_max; ++pat_one_bt_cnt) {
1888 Segment = (EdbSegP*)pat_one->GetSegment(pat_one_bt_cnt);
1889
1890 IP_Pair_To_InBT_Seg =CalcIP(Segment,vtx);
1891 IP_Pair_To_InBT=IP_Pair_To_InBT_Seg;
1892 IP_Seg1_To_Vtx = CalcIP(Segment,vtx);
1893
1894 // Check if IP to vtx (BT) is ok:
1895 if (IP_Pair_To_InBT>eParaValue[0]) continue;
1896
1897 //-----------------------------------
1898 // --- Loop over basetracks of second pattern
1899 for (Int_t pat_two_bt_cnt=0; pat_two_bt_cnt<pat_two_bt_cnt_max; ++pat_two_bt_cnt) {
1900 Segment2 = (EdbSegP*)pat_two->GetSegment(pat_two_bt_cnt);
1901
1902 // For statistics: increase the number of total searched combinations:
1903 ++NPairTriesTotal;
1904
1905 if (Segment2==Segment) continue;
1906 if (Segment2->ID()==Segment->ID()&&Segment2->PID()==Segment->PID()) continue;
1907
1908 // Now apply cut conditions: GS GAMMA SEARCH Alg --------------------
1909 // if InBT is flagged as MC InBT, take care that only BG or same MC basetracks are taken:
1910 if (InBT->MCEvt()>0) if (Segment->MCEvt()>0&&Segment2->MCEvt()>0) if (Segment->MCEvt()!=Segment2->MCEvt()) continue;
1911 if (InBT->MCEvt()>0) if (Segment->MCEvt()>0&&Segment2->MCEvt()>0) if (Segment->MCEvt()!=InBT->MCEvt()) continue;
1912 if (InBT->MCEvt()>0) if (Segment->MCEvt()>0&&Segment2->MCEvt()>0) if (Segment2->MCEvt()!=InBT->MCEvt()) continue;
1913
1914 // In case of two MC events, check for e+ e- pairs
1915 // Only if Parameter is set to choose different Flag() pairs:
1916 if (InBT->MCEvt()>0 && eParaValue[6]==1)
1917 if (Segment->MCEvt()>0&&Segment2->MCEvt()>0)
1918 if ((Segment2->Flag()+Segment->Flag())!=0) continue;
1919
1920 // a) Check dR between tracks:
1921 if (DeltaR_WithPropagation(Segment,Segment2)>eParaValue[2]) continue;
1922 // b) Check dT between tracks:
1923 if (DeltaThetaSingleAngles(Segment,Segment2)>eParaValue[4]) continue;
1924 // c) Check dMinDist between tracks:
1925 if (GetMinimumDist(Segment,Segment2)>eParaValue[1]) continue;
1926
1927 // d) Check if dist Z to vtx (BT) is ok:
1928 distZ=Segment->Z()-InBT->Z();
1929 if (distZ>eParaValue[3]) continue;
1930
1931 // At first: Check for already duplicated pairings:
1932 if (CheckPairDuplications(Segment->PID(),Segment->ID(),Segment2->PID(),Segment2->ID(), eSegmentPIDArray,eSegmentIDArray,eSegment2PIDArray,eSegment2IDArray,RecoShowerArrayN)) continue;
1933
1934 // Check if both basetracks have a vertex which is upstream
1935 // of both tracks (only then the two BT are really pointing).
1936 TObjArray *segments = new TObjArray(2);
1937 segments->Add(Segment);
1938 segments->Add(Segment2);
1939 EdbVertex* vetex = CalcVertex(segments);
1940 if (vetex ->Z()> TMath::Min(Segment->Z(),Segment2->Z()) ) continue;
1941
1942 // Set Sum values of the to BTs:
1943 x_av=Segment2->X()+(Segment->X()-Segment2->X())/2.0;
1944 y_av=Segment2->Y()+(Segment->Y()-Segment2->Y())/2.0;
1945 z_av=Segment2->Z()+(Segment->Z()-Segment2->Z())/2.0;
1946 tx_av=Segment2->TX()+(Segment->TX()-Segment2->TX())/2.0;
1947 ty_av=Segment2->TY()+(Segment->TY()-Segment2->TY())/2.0;
1948 Segment_Sum->SetX(x_av);
1949 Segment_Sum->SetY(y_av);
1950 Segment_Sum->SetTX(tx_av);
1951 Segment_Sum->SetTY(ty_av);
1952 Segment_Sum->SetZ(z_av);
1953
1954 // Check if IP to vtx (BT) is ok:
1955 IP_Pair_To_InBT_Seg2 = CalcIP(Segment2,vtx);
1956 if (IP_Pair_To_InBT_Seg2>eParaValue[0]) continue;
1957
1958 IP_Seg2_To_Vtx = CalcIP(Segment2,vtx);
1959 IP_Pair_To_InBT_SegSum=CalcIP(Segment_Sum, vtx);
1960
1961 Float_t IP_Seg1ToVtxSeg1Seg2=0;
1962 Float_t IP_Seg2ToVtxSeg1Seg2=0;
1963 IP_Seg1ToVtxSeg1Seg2 = CalcIP(Segment ,vetex);
1964 IP_Seg2ToVtxSeg1Seg2 = CalcIP(Segment2,vetex);
1965
1966 if (gEDBDEBUGLEVEL>2) {
1967 cout << "EdbShowerAlg_GS::FindPairs IP_Pair_To_InBT_Seg = " << IP_Pair_To_InBT_Seg << endl;
1968 cout << "EdbShowerAlg_GS::FindPairs IP_Pair_To_InBT_Seg2 = " << IP_Pair_To_InBT_Seg2 << endl;
1969 cout << "EdbShowerAlg_GS::FindPairs IP_Pair_To_InBT_SegSum= " << IP_Pair_To_InBT_SegSum << endl;
1970 cout << "EdbShowerAlg_GS::FindPairs IP_Seg1ToVtxSeg1Seg2= " << IP_Seg1ToVtxSeg1Seg2 << endl;
1971 cout << "EdbShowerAlg_GS::FindPairs IP_Seg2ToVtxSeg1Seg2= " << IP_Seg2ToVtxSeg1Seg2 << endl;
1972 }
1973
1974 // Save the segment which has smaller IP, this will be the first BT in the RecoShower
1975 if ( IP_Pair_To_InBT_Seg>IP_Pair_To_InBT_Seg2 ) {
1976 IP_Pair_To_InBT_SegSmaller=0; // take Segment first
1977 }
1978 else {
1979 IP_Pair_To_InBT_SegSmaller=1; // take Segment2 first
1980 }
1981
1982 // e) Check if this is not a possible fake doublet which is
1983 // sometimes caused by view overlap in the scanning:
1984 // Deprecated, because this is now incorporated in the cleaning routine.
1985 //
1986 // f) new: Last check if ePair is from e+e-Pair or from fake BG;
1987 // do a NN calculation estimate:
1988 //cout << "DEBUGTEST: PRINTNICE SEGMENT AND SEGMENT2 AND VTX" << endl;
1989 //Segment->PrintNice();
1990 //Segment2->PrintNice();
1991 //cout << "Vtx: Print MC XYZ: " << vtx ->MCEvt() << " " << vtx ->X() << " " << vtx ->Y() << " " << vtx ->Z() <<endl;
1992
1993 eValueGSNN_varInput=-1;
1994 // Purity 1: Input = 1.0;
1995 // Purity 0.5: Input = 0.5;
1996 // Purity else: Input = 0.0;
1997 if (Segment2->Flag()+Segment->Flag()==0&&TMath::Abs(Segment2->Flag())==11&&Segment->MCEvt()>0) {
1998 eValueGSNN_varInput=1;
1999 }
2000 else if (Segment2->Flag()+Segment->Flag()!=0&&TMath::Abs(Segment2->Flag())==11&&Segment->MCEvt()>0) {
2001 eValueGSNN_varInput=0.5;
2002 }
2003 else {
2004 eValueGSNN_varInput=0;
2005 }
2006
2007 eValueGSNN_var00=TMath::Min(IP_Seg1_To_Vtx,IP_Seg2_To_Vtx);
2008 if (eRecoMode==1) eValueGSNN_var00=IP_InBT_To_Vtx;
2009 if (eRecoMode==2) {
2010 eValueGSNN_var00=TMath::Min(IP_Seg1ToVtxSeg1Seg2,IP_Seg2ToVtxSeg1Seg2);
2011 }
2012
2013 // Set Tree Input Variable Values:
2014 eValueGSNN_var01=GetMinimumDist(Segment,Segment2);
2015 eValueGSNN_var02=DeltaR_WithPropagation(Segment,Segment2);
2016 eValueGSNN_var03=InBT->Z()-vtx->Z();
2017 if (eRecoMode==2) {
2018 eValueGSNN_var03=TMath::Min(Segment->Z(),Segment2->Z())-vetex->Z();
2019 }
2020 if (eRecoMode==0) {
2021 eValueGSNN_var03=TMath::Min(Segment->Z(),Segment2->Z())-vtx->Z();
2022 }
2023 eValueGSNN_var04=DeltaThetaSingleAngles(Segment,Segment2);
2024 eValueGSNN_var05=TMath::Abs(pat_one_cnt-pat_two_cnt);
2025 eValueGSNN_var06=Segment2->Flag()+Segment->Flag();
2026
2027 if (gEDBDEBUGLEVEL>2) {
2028 cout <<"------------- Print eValueGSNN ------------- " << endl;
2029 Segment->PrintNice();
2030 Segment2->PrintNice();
2031 cout <<"------------- Print eValueGSNN_var00= " << eValueGSNN_var00 << " ---- " << endl;
2032 cout <<"------------- Print eValueGSNN_var01= " << eValueGSNN_var01 << " ---- " << endl;
2033 cout <<"------------- Print eValueGSNN_var02= " << eValueGSNN_var02 << " ---- " << endl;
2034 cout <<"------------- Print eValueGSNN_var03= " << eValueGSNN_var03 << " ---- " << endl;
2035 cout <<"------------- Print eValueGSNN_var04= " << eValueGSNN_var04 << " ---- " << endl;
2036 cout <<"------------- Print eValueGSNN_var05= " << eValueGSNN_var05 << " ---- " << endl;
2037 }
2038
2039 // Severe Warning! If in this part of the code:
2040 // eANNPairTree->Show(eANNPairTree->GetEntries()-1) is invoked,
2041 // then the Variables of the Treeentry are loaded into the memory
2042 // address and the current value is overwritten!!!!!!
2043 // DO NOT DO THIS!
2044
2045 // Now we have to fill params values.
2046 // Distinguish CASE A,C and CASE B, because of one less input neuron.
2047 if (eRecoMode==0||eRecoMode==2) {
2048 paramsAC[0]=eValueGSNN_var00; //dIP
2049 paramsAC[1]=eValueGSNN_var01; //dMinDist
2050 paramsAC[2]=eValueGSNN_var02; //dR
2051 paramsAC[3]=eValueGSNN_var03; //deltaZ
2052 paramsAC[4]=eValueGSNN_var04; //dT
2053 paramsAC[5]=eValueGSNN_var05; //dNPL
2054 //paramsAC[6]=eValueGSNN_var06; // cannot be used, uses MC info...
2055 }
2056 else {
2057 paramsB[0]=eValueGSNN_var00; //dIP
2058 paramsB[0]=eValueGSNN_var01; //dMinDist
2059 paramsB[1]=eValueGSNN_var02; //dR
2060 //paramsB[3]=eValueGSNN_var03;
2061 //deltaZ // constant, so take out (otherwise ANN doesnt work!)
2062 paramsB[2]=eValueGSNN_var04; //dT
2063 paramsB[3]=eValueGSNN_var05; //dNPL
2064 }
2065
2066 //---------
2067 Double_t value=0;
2068 if (gEDBDEBUGLEVEL>2) {
2069 int nparams=6;
2070 if (eRecoMode==1) {
2071 nparams=4;
2072 cout << "EdbShowerAlg_GS::FindPairs ParamsB: ";
2073 for (int hj=0; hj<nparams; hj++) cout << " " << paramsB[hj];
2074 cout << " - Input: " << eValueGSNN_varInput << endl;
2075 cout << "EdbShowerAlg_GS::FindPairs value=eANNPair->Evaluate(0,params) = ";
2076 value=eANNPair->Evaluate(0,paramsB);
2077 eValueGSNN_varOutput=value;
2078 cout << value << endl;
2079 cout << "EdbShowerAlg_GS::FindPairs -------- " << endl;
2080 }
2081 else {
2082 nparams=6;
2083 cout << "EdbShowerAlg_GS::FindPairs ParamsAC: ";
2084 for (int hj=0; hj<nparams; hj++) cout << " " << paramsAC[hj];
2085 cout << " - Input: " << eValueGSNN_varInput << endl;
2086 cout << "EdbShowerAlg_GS::FindPairs value=eANNPair->Evaluate(0,params) = ";
2087 value=eANNPair->Evaluate(0,paramsAC);
2088 eValueGSNN_varOutput=value;
2089 cout << value << endl;
2090 cout << "EdbShowerAlg_GS::FindPairs -------- " << endl;
2091 }
2092 }
2093
2094
2095 // Evaluate now the NN.
2096 if (eRecoMode==1) {
2097 value=eANNPair->Evaluate(0,paramsB);
2098 }
2099 else {
2100 value=eANNPair->Evaluate(0,paramsAC);
2101 }
2102
2103
2104 eValueGSNN_varOutput=value;
2105 // Fill the Tree with NN variables later. Then it can be done
2106 // also for the cleaned pairs.
2107 // eANNPairTree->Fill();
2108
2109 // After Full Tree is Filled, check wheter last cut condition is satisfied.
2110 if ( eCutModeFull ) if (value<eANNPairCut[eRecoMode]) continue;
2111 if (value<eANNPairCut[eRecoMode]) continue;
2112
2113
2114// if (gEDBDEBUGLEVEL>2) {
2115 if (value<eANNPairCut[eRecoMode]) {
2116 cout << "EdbShowerAlg_GS::FindPairs This Pair has passed all cuts w.r.t to InBT:" << endl;
2117 cout << "EdbShowerAlg_GS::FindPairs IP_Pair_To_InBT = " << IP_Pair_To_InBT << endl;
2118 cout << "EdbShowerAlg_GS::FindPairs GetMinimumDist(Segment,Segment2) = " << GetMinimumDist(Segment,Segment2) << endl;
2119 cout << "EdbShowerAlg_GS::FindPairs CalcIP(Segment_Sum,InBT) = " << IP_Pair_To_InBT << endl;
2120 cout << "EdbShowerAlg_GS::FindPairs IP_Pair_To_InBT_Seg = " << IP_Pair_To_InBT_Seg << endl;
2121 cout << "EdbShowerAlg_GS::FindPairs IP_Pair_To_InBT_Seg2 = " << IP_Pair_To_InBT_Seg2 << endl;
2122 cout << "EdbShowerAlg_GS::FindPairs GetSpatialDist(Segment_Sum,InBT) = " << GetSpatialDist(Segment_Sum,InBT) << endl;
2123 cout << "EdbShowerAlg_GS::FindPairs IP_InBT_To_Vtx = " << IP_InBT_To_Vtx << endl;
2124 cout << "EdbShowerAlg_GS::FindPairs IP_Seg1_To_Vtx = " << IP_Seg1_To_Vtx << endl;
2125 cout << "EdbShowerAlg_GS::FindPairs IP_Seg2_To_Vtx = " << IP_Seg2_To_Vtx << endl;
2126 cout << "InBT: Print InBT: " << InBT ->X() << " " << InBT ->Y() << " " << InBT ->Z() <<endl;
2127 cout << "Vtx: Print MC XYZ: " << vtx ->MCEvt() << " " << vtx ->X() << " " << vtx ->Y() << " " << vtx ->Z() <<endl;
2128 }
2129
2130 eSegmentPIDArray->AddAt(Segment->PID(), RecoShowerArrayN);
2131 eSegmentIDArray->AddAt(Segment->ID(), RecoShowerArrayN);
2132 eSegment2PIDArray->AddAt(Segment2->PID(), RecoShowerArrayN);
2133 eSegment2IDArray->AddAt(Segment2->ID(), RecoShowerArrayN);
2134
2135
2136 // Create new EdbTrackP Object for storage;
2137 EdbTrackP* RecoShower = new EdbTrackP();
2138 RecoShower -> AddSegment(Segment);
2139 RecoShower -> AddSegment(Segment2);
2140 // Set X and Y and Z values: /Take lower Z of both BTs)
2141 RecoShower->SetZ(TMath::Min(Segment->Z(),Segment2->Z()));
2142 RecoShower->SetX(Segment_Sum->X());
2143 RecoShower->SetY(Segment_Sum->Y());
2144 RecoShower->SetTX(Segment_Sum->TX());
2145 RecoShower->SetTY(Segment_Sum->TY());
2146 RecoShower->SetMC(InBT->MCEvt(),InBT->MCEvt());
2147 RecoShower->SetID(eRecoShowerArrayN);
2148 RecoShower->SetPID(Segment->PID());
2149
2150 // Add Shower to interim Array:
2151 RecoShowerArray->Add(RecoShower);
2152 ++RecoShowerArrayN;
2153
2154 if (gEDBDEBUGLEVEL>2) {
2155 cout << "Added shower at " << RecoShower << " to reco shower array. Until now: ";
2156 cout << RecoShowerArrayN << " entries in it: Print Shower:" << endl;
2157 RecoShower->PrintNice();
2158 }
2159
2160 delete segments;
2161 }
2162 }
2163 }
2164 } //for (Int_t pat_one_cnt=0; ...
2165
2166 // ------------------------------------
2167 if (gEDBDEBUGLEVEL>2) {
2168 cout << "EdbShowerAlg_GS::FindPairs For the InBT/Vtx at __" << InBT << "__, we have searched " << NPairTriesTotal << " pair combinations in the PVRec volume." << endl;
2169 cout << "EdbShowerAlg_GS::FindPairs For the InBT/Vtx at __" << InBT << "__, we have found " << RecoShowerArray->GetEntries() << " compatible pairs in the PVRec volume." << endl;
2170 }
2171 // ------------------------------------
2172
2173 // Delete unnecessary objects: important, else memory overflow!
2174 delete vtx;
2175 delete Segment_Sum;
2176
2177 Log(3,"EdbShowerAlg_GS::FindPairs","Starting FindPairs(InBT,eAli_Sub) now....done.");
2178 return RecoShowerArray;
2179}
npat
Int_t npat
Definition: Xi2HatStartScript.C:33
EdbPattern
Definition: EdbPattern.h:273
EdbPatternsVolume::GetPatternZLowestHighest
EdbPattern * GetPatternZLowestHighest(Bool_t lowestZ=kTRUE) const
Definition: EdbPattern.cxx:1816
EdbPatternsVolume::Npatterns
Int_t Npatterns() const
Definition: EdbPattern.h:366
EdbPatternsVolume::GetPattern
EdbPattern * GetPattern(int id) const
Definition: EdbPattern.cxx:1721
EdbSegP::SetPID
void SetPID(int pid)
Definition: EdbSegP.h:129
EdbSegP::SetY
void SetY(Float_t y)
Definition: EdbSegP.h:178
EdbSegP::TX
Float_t TX() const
tangens = deltaX/deltaZ
Definition: EdbSegP.h:175
EdbSegP::SetTX
void SetTX(Float_t tx)
Definition: EdbSegP.h:179
EdbSegP::ID
Int_t ID() const
Definition: EdbSegP.h:147
EdbSegP::SetID
void SetID(int id)
Definition: EdbSegP.h:128
EdbSegP::SetX
void SetX(Float_t x)
Definition: EdbSegP.h:177
EdbSegP::X
Float_t X() const
Definition: EdbSegP.h:173
EdbSegP::SetTY
void SetTY(Float_t ty)
other functions
Definition: EdbSegP.h:180
EdbSegP::Y
Float_t Y() const
Definition: EdbSegP.h:174
EdbSegP::TY
Float_t TY() const
tangens = deltaY/deltaZ
Definition: EdbSegP.h:176
EdbSegmentsBox::Z
Float_t Z() const
Definition: EdbPattern.h:84
EdbSegmentsBox::Y
Float_t Y() const
Definition: EdbPattern.h:83
EdbSegmentsBox::GetN
Int_t GetN() const
Definition: EdbPattern.h:65
EdbSegmentsBox::GetSegment
EdbSegP * GetSegment(int i) const
Definition: EdbPattern.h:66
EdbSegmentsBox::X
Float_t X() const
Definition: EdbPattern.h:82
EdbShowerAlg_GS::CheckPairDuplications
Bool_t CheckPairDuplications(Int_t SegPID, Int_t SegID, Int_t Seg2PID, Int_t Seg2ID, TArrayI *SegmentPIDArray, TArrayI *SegmentIDArray, TArrayI *Segment2PIDArray, TArrayI *Segment2IDArray, Int_t RecoShowerArrayN)
Definition: EdbShowerAlg.cxx:1158
EdbShowerAlg_GS::eANNPairCut
Float_t eANNPairCut[3]
Definition: EdbShowerAlg.h:221
EdbShowerAlg_GS::eCutModeFull
Bool_t eCutModeFull
Definition: EdbShowerAlg.h:206
EdbShowerAlg::GetSpatialDist
Double_t GetSpatialDist(EdbSegP *s1, EdbSegP *s2)
Definition: EdbShowerAlg.cxx:633
EdbShowerAlg::eParaValue
Float_t eParaValue[10]
Definition: EdbShowerAlg.h:41
MIPEmulsionTrack::Z
Double_t Z
Definition: tlg2couples.C:104

◆ FindPairsPreselected() [1/2]

TObjArray * EdbShowerAlg_GS::FindPairsPreselected ( EdbPVRec *  eAli_Sub)

FindPairsPreselected()
What is the goal of this function?
A simple example:
Consider 30 plates, Each plate has 20000 basetracks in it.
GammaPairing from volume would then have to search all these
30*20000=600000 BTs for combinations.
That means that for each of these 600k basetracks, we would have
to loop and compare with 4*20000 other basetracks to see if we
have BT pair combinations possible this each basetrack.
(4 is the maximum difference number in plates according to this algo).
So 80000 possible combinations to check for one BT, makes in total

80000*600000 bt pairs. That is too much for reasonable time. \n

Small Trick.... The BTs are sorted w.r.t. z position!

1526{
1540
1541 Log(2,"EdbShowerAlg_GS::FindPairsPreselected","FindPairsPreselected()...");
1542
1543 EdbSegP* InBT=NULL;
1544 EdbSegP* Segment=NULL;
1545 EdbSegP* Segment2=NULL;
1546 TObjArray* PreselectedArray= new TObjArray();
1547 TObjArray* NewPreselectedArray= new TObjArray();
1548 Int_t counter=0;
1549 Int_t newcounter=0;
1550 Int_t npat=eAli_Sub->Npatterns();
1551 Int_t pat_one_bt_cnt_max,pat_two_bt_cnt_max=0;
1552 EdbPattern* pat_one=0;
1553 EdbPattern* pat_two=0;
1554 EdbPVRec* pvr=eAli_Sub;
1555
1556 cout << "EdbShowerAlg_GS::FindPairsPreselected eAli_Sub->GetEntries() " << eAli_Sub->NSeg() << endl;
1558
1559 for (int l=0; l<eAli_Sub->Npatterns(); l++) {
1560 if (gEDBDEBUGLEVEL>2) cout << "Doing pattern() " << l << endl;
1561
1562 // EdbPattern with couples
1563 EdbPattern *pat = pvr->GetPattern(l);
1564
1565 // Make a hash table in X,Theta space, with division of x microns and y rad.
1566 pat->FillCell(500, 500, 0.1, 0.1);
1567
1568 // Loop over all segments.
1569 for (int j=0; j<pat->N(); j++) {
1570 EdbSegP *s1 = pat->GetSegment(j);
1571 // s1 is used as prediction. give your criteria as one sigma.
1572 // set value should be square. (dx^2, dy^2, dz^2, dtx^2, dty^2)
1573 s1->SetErrors(3000,3000,100, 0.05,0.05);
1574 // dx=10microns, dy=10microns, dz=sqrt(10)microns, dtx=0.1rad, dty=0.1rad.
1575 // get corresponding tracks in pattern.
1576 // pat->FindCompliments() will search segments in hash table and put them in segments.
1577 // last 2 numbers are number of sigma to be taken.
1578 TObjArray segments; // array to be filled with segments.
1579 Int_t ncompl = pat->FindCompliments(*s1, segments, 1, 1);
1580 //cout << "j= " << j << " ncompl= " << ncompl << " segments.GetEntries() " << segments.GetEntries() << endl;
1581
1582 // ncompl==1 means only the BT itself was found:
1583 if (ncompl==1) continue;
1584
1585 for (int k=0; k<segments.GetEntries(); k++) {
1586 EdbSegP *s2 = (EdbSegP *) segments.At(k);
1587 NewPreselectedArray->Add(s2);
1588 }
1589 }
1590 // for (int j=0; j<pat->N(); j++)
1591 } // for (int l=0; l<eAli_Sub->Npatterns(); l++)
1592 cout << "EdbShowerAlg_GS::FindPairsPreselected Preselected: "<< NewPreselectedArray->GetEntries() << endl;
1593
1594 // Here it might be that there are still duplicated basetracks (per plate) in the array,
1595 // these we have to take out:
1596 cout << "EdbShowerAlg_GS::FindPairsPreselected Check for duplicated basetracks:" << endl;
1597 cout << "EdbShowerAlg_GS::FindPairsPreselected This might take a while...:" << endl;
1598 cout << "EdbShowerAlg_GS::FindPairsPreselected (status bar):" << endl;
1599 Int_t count=0;
1600 Int_t Zseparator[114];
1601 for (Int_t l=0; l< 114; l++) {
1602 Zseparator[l]=0;
1603 }
1604 EdbSegP *sl = (EdbSegP*) NewPreselectedArray->At(0);
1605 Zseparator[0]=0;
1606 Float_t ZStart=sl->Z();
1607 Int_t incrementor=1;
1608 Int_t NewPreselectedArrayN=NewPreselectedArray->GetEntries();
1609
1611 for (Int_t l=0; l< NewPreselectedArrayN; l++) {
1612 EdbSegP *sl = (EdbSegP*) NewPreselectedArray->At(l);
1613 if (sl->Z()!=ZStart) {
1614 Zseparator[incrementor]=l;
1615 ZStart=sl->Z();
1616 incrementor++;
1617 }
1618 }
1619
1620
1621 if (gEDBDEBUGLEVEL>2) {
1622 for (Int_t l=0; l< 114; l++) {
1623 cout << "EdbShowerAlg_GS::FindPairsPreselected incrementor = " << incrementor << endl;
1624 }
1625 }
1626
1627
1628 for (Int_t o=0; o < incrementor-1; o++) {
1629
1630 cout << "." ;
1631
1632 if (gEDBDEBUGLEVEL>2) cout << "EdbShowerAlg_GS::FindPairsPreselected Doing ZSeparator o= " << o << " loop from BT #" << Zseparator[o] << " to BT # " << Zseparator[o+1]-1 << endl;
1633
1634
1635 for (Int_t l=Zseparator[o]; l< Zseparator[o+1]-1; l++) {
1636
1637 Bool_t IsDuplicateThere=kFALSE;
1638 EdbSegP *sl = (EdbSegP*) NewPreselectedArray->At(l);
1639
1640
1641 for (Int_t k=l+1; k< Zseparator[o+1]; k++) {
1642 EdbSegP *sk = (EdbSegP*) NewPreselectedArray->At(k);
1643
1644 if (TMath::Abs(sl->X()-sk->X())<0.1 && TMath::Abs(sl->Y()-sk->Y())<0.1 ) {
1645 if (gEDBDEBUGLEVEL>2) {
1646 cout << "Duplicate found: " << sl << " and " << sk << endl;
1647 sl->PrintNice();
1648 sk->PrintNice();
1649 }
1650 ++count;
1651 IsDuplicateThere=kTRUE;
1652 }
1653 } //for (Int_t k=l+1; k< Zseparator[o+1]; k++)
1654 if (!IsDuplicateThere) PreselectedArray->Add(sl);
1655
1656 } // for (Int_t l=Zseparator[o]; l< Zseparator[o+1]-1; l++)
1657 } // for (Int_t o=0; o < incrementor; o++) {
1658
1659
1660 if (gEDBDEBUGLEVEL>2) {
1661 cout << "EdbShowerAlg_GS::FindPairsPreselected NDuplicates = " << count << endl;
1662 cout << "EdbShowerAlg_GS::FindPairsPreselected NNewPreselectedArray = " << NewPreselectedArray->GetEntries() << endl;
1663 }
1664 cout << "EdbShowerAlg_GS::FindPairsPreselected NPreselectedArray = " << PreselectedArray->GetEntries() << endl;
1665
1666 Log(2,"EdbShowerAlg_GS::FindPairsPreselected","FindPairsPreselected()...done.");
1667 return PreselectedArray;
1669}
EdbPVRec
Definition: EdbPVRec.h:148
EdbPVRec::NSeg
Int_t NSeg()
Definition: EdbPVRec.cxx:2913
EdbPattern::FindCompliments
int FindCompliments(EdbSegP &s, TObjArray &arr, float nsig, float nsigt)
Definition: EdbPattern.cxx:1354
EdbPattern::FillCell
void FillCell(float stepx, float stepy, float steptx, float stepty)
Definition: EdbPattern.cxx:1323
EdbSegP::SetErrors
void SetErrors()
Definition: EdbSegP.h:90
EdbSegmentsBox::N
Int_t N() const
Definition: EdbPattern.h:86
o
AcqOdyssey * o
Definition: hwinit.C:2

◆ FindPairsPreselected() [2/2]

Bool_t EdbShowerAlg_GS::FindPairsPreselected ( EdbSegP *  InitiatorBT,
EdbPVRec *  eAli_Sub 
)
1675{
1676 Log(3,"EdbShowerAlg_GS::FindPairsPreselected","Starting FindPairsPreselected(InBT,eAli_Sub) now....");
1677
1678 EdbSegP* InBT=NULL;
1679 InBT=InitiatorBT;
1680 EdbSegP* Segment=NULL;
1681 EdbSegP* Segment2=NULL;
1682
1683 Int_t npat=eAli_Sub->Npatterns();
1684 Int_t pat_one_bt_cnt_max,pat_two_bt_cnt_max=0;
1685 EdbPattern* pat_one=0;
1686 EdbPattern* pat_two=0;
1687
1688 //-----------------------------------
1689 // --- Loop over first patterns
1690 for (Int_t pat_one_cnt=0; pat_one_cnt<npat; ++pat_one_cnt) {
1691 pat_one=(EdbPattern*)eAli_Sub->GetPattern(pat_one_cnt);
1692 pat_one_bt_cnt_max=eAli_Sub->GetPattern(pat_one_cnt)->GetN();
1693
1694 // Check if dist Z to vtx (BT) is ok:
1695 Float_t distZ=pat_one->Z()-InBT->Z();
1696 if (gEDBDEBUGLEVEL>4) {
1697 cout << "EdbShowerAlg_GS::FindPairsPreselected Check if dist Z to vtx (BT) is ok: distZ=" << distZ << endl;
1698 }
1699 if (distZ>eParaValue[3]) continue;
1700 if (distZ<-1000) continue;
1701
1702 if (gEDBDEBUGLEVEL>2) cout << "EdbShowerAlg_GS::FindPairsPreselected pat_one_cnt=" << pat_one_cnt << " pat_one->Z() = " << pat_one->Z() << " pat_one_bt_cnt_max= "<< pat_one_bt_cnt_max <<endl;
1703
1704 //-----------------------------------
1705 // --- Loop over second patterns
1706 for (Int_t pat_two_cnt=0; pat_two_cnt<npat; ++pat_two_cnt) {
1707
1708 // Now apply cut conditions: GS GAMMA SEARCH Alg:
1709 if (TMath::Abs(pat_one_cnt-pat_two_cnt)>eParaValue[5]) continue;
1710
1711 pat_two=(EdbPattern*)eAli_Sub->GetPattern(pat_two_cnt);
1712 pat_two_bt_cnt_max=eAli_Sub->GetPattern(pat_two_cnt)->GetN();
1713
1714 distZ=pat_two->Z()-InBT->Z();
1715 if (distZ>eParaValue[3]) continue;
1716 if (distZ<-1000) continue;
1717
1718 if (gEDBDEBUGLEVEL>2) cout << "EdbShowerAlg_GS::FindPairsPreselected pat_two_cnt=" << pat_two_cnt << " pat_two->Z() = " << pat_two->Z() << " pat_two_bt_cnt_max= "<< pat_two_bt_cnt_max <<endl;
1719
1720 //-----------------------------------
1721 // --- Loop over basetracks of first pattern
1722 for (Int_t pat_one_bt_cnt=0; pat_one_bt_cnt<pat_one_bt_cnt_max; ++pat_one_bt_cnt) {
1723 Segment = (EdbSegP*)pat_one->GetSegment(pat_one_bt_cnt);
1724
1725 //-----------------------------------
1726 // --- Loop over basetracks of second pattern
1727 for (Int_t pat_two_bt_cnt=0; pat_two_bt_cnt<pat_two_bt_cnt_max; ++pat_two_bt_cnt) {
1728
1729 Segment2 = (EdbSegP*)pat_two->GetSegment(pat_two_bt_cnt);
1730 if (Segment2==Segment) continue;
1731 if (Segment2->ID()==Segment->ID()&&Segment2->PID()==Segment->PID()) continue;
1732
1733 // a) Check dR between tracks:
1734 if (DeltaR_WithPropagation(Segment,Segment2)<eParaValue[2]) return kTRUE;
1735 // b) Check dT between tracks:
1736 if (DeltaThetaSingleAngles(Segment,Segment2)>eParaValue[4]) return kTRUE;
1737 // c) Check dMinDist between tracks:
1738 if (GetMinimumDist(Segment,Segment2)>eParaValue[1]) return kTRUE;
1739
1740 }
1741 }
1742 }
1743 } //for (Int_t pat_one_cnt=0; ...
1744
1745 Log(3,"EdbShowerAlg_GS::FindPairsPreselected","Starting FindPairsPreselected(InBT,eAli_Sub) now....done.");
1746 return kFALSE;
1747}

◆ GetInVtxArray()

TObjArray * EdbShowerAlg_GS::GetInVtxArray ( ) const
inline
259 {
260 return eInVtxArray;
261 }

◆ GetInVtxArrayN()

Int_t EdbShowerAlg_GS::GetInVtxArrayN ( ) const
inline
256 {
257 return eInVtxArrayN;
258 }

◆ Init()

void EdbShowerAlg_GS::Init ( void  )

Init with values according to GS Alg.
Create essential objects.
Init function is supposed to be executed exact once for
each class instance.
Took over from "FindGamma.C" script I develoved before:
and with recorded values from a GammaSearch best algo parameterset.
See also table "tableshowerchosenparametersets"
and figure "" in FWM thesis.

812{
823
824 Log(2,"EdbShowerAlg_GS::EdbShowerAlg_GS","Init()");
825
826 eDebug=kFALSE;
827
828 eRecoMode=0;
829 eCutModeFull=kTRUE;
830
831 eAlgName="GS";
832 eAlgValue=999; // see default.par_SHOWREC for labeling (labeling identical with ShowRec program)
833
834 eFindPairsPreselected=kTRUE;
835
836 eANNPairCut[0]=0.5;
837 eANNPairCut[1]=0.48;
838 eANNPairCut[2]=0.45;
839
840 //
841 // min IP cut; this cut is used for the -better- IP of both BTs to Vertex/BT
842 eParaValue[0]=393;
843 eParaString[0]="PARA_GS_CUT_dIP";
844 //
845 // min minDist.e between pair BTs
846 eParaValue[1]=35;
847 eParaString[1]="PARA_GS_CUT_dMin";
848 //
849 // min dR (with propagation) between pair BTs
850 eParaValue[2]=85;
851 eParaString[2]="PARA_GS_CUT_dR";
852 //
853 // max Z distance between pair BTs and Vertex/BT
854 eParaValue[3]=25000;
855 eParaString[3]="PARA_GS_CUT_dZ";
856 //
857 // max Angle between pair BTs
858 eParaValue[4]=0.11;
859 eParaString[4]="PARA_GS_CUT_dtheta";
860 //
861 // max plates difference between pair BTs
862 eParaValue[5]=1;
863 eParaString[5]="PARA_GS_CUT_PIDDIFF";
864 //
865 // in MC case: check for opposite flag sign
866 eParaValue[6]=1;
867 eParaString[6]="PARA_GS_CUT_OPPOSITEFLAG";
868
869
870 // Variable to clean the found parings once more.
871 // Default is kTRUE.
872 eSetCleanPairs=kFALSE;
873 eSetCleanPairs=kTRUE;
874
875 // Create needed objects for storing Initiator Vertices/Basetracks.
876 eRecoShowerArray= new TObjArray(999);
877 eRecoShowerArrayN=0;
878
879 eInVtxArray = new TObjArray();
880 eInVtxArrayN=0;
881
882 eInBTArray = new TObjArray();
883 eInBTArrayN=0;
884
885 eSegmentPIDArray = new TArrayI(9999);
886 eSegmentIDArray = new TArrayI(9999);
887 eSegment2PIDArray = new TArrayI(9999);
888 eSegment2IDArray = new TArrayI(9999);
889
890 if (gEDBDEBUGLEVEL>2) {
891 cout << "EdbShowerAlg_GS eRecoShowerArray = " << eRecoShowerArray << endl;
892 cout << "EdbShowerAlg_GS eInVtxArray = " << eInVtxArray << endl;
893 cout << "EdbShowerAlg_GS eInBTArray = " << eInBTArray << endl;
894 }
895
896 Log(2,"EdbShowerAlg_GS::EdbShowerAlg_GS","Init()...done.");
897 return;
898}
EdbShowerAlg::eParaString
TString eParaString[10]
Definition: EdbShowerAlg.h:42
EdbShowerAlg::eAlgValue
Int_t eAlgValue
Definition: EdbShowerAlg.h:40
EdbShowerAlg::eAlgName
TString eAlgName
Definition: EdbShowerAlg.h:39

◆ Initialize()

void EdbShowerAlg_GS::Initialize ( )
virtual

Reimplemented from EdbShowerAlg.

1024{
1025 Log(2,"EdbShowerAlg_GS::EdbShowerAlg_GS","Initialize()");
1026 return;
1027}

◆ IsPossibleFakeDoublet()

Bool_t EdbShowerAlg_GS::IsPossibleFakeDoublet ( EdbSegP *  s1,
EdbSegP *  s2 
)
2422{
2423 if (TMath::Abs(s1->X()-s2->X())<1) return kTRUE; // minimum distance of 1micron
2424 if (TMath::Abs(s1->Y()-s2->Y())<1) return kTRUE;// minimum distance of 1micron
2425 if (TMath::Abs(s1->TX()-s2->TX())<0.005) return kTRUE;// minimum angle of 5mrad
2426 if (TMath::Abs(s1->TY()-s2->TY())<0.005) return kTRUE;// minimum angle of 5mrad
2427 return kFALSE;
2428}

◆ ReloadANNs()

void EdbShowerAlg_GS::ReloadANNs ( Int_t  RecoMode)
994{
995 Log(2,"EdbShowerAlg_GS::ReloadANNs","ReloadANNs(Int_t RecoMode)");
996
997 if (RecoMode==0) {
998 eANNPair=eANNPairCaseA;
999 Log(2,"EdbShowerAlg_GS::ReloadANNs","Set eANNPair to eANNPairCaseA");
1000 cout << "EdbShowerAlg_GS::ReloadANNs: eANNPairCaseA at " << eANNPairCaseA << endl;
1001 }
1002 else if (RecoMode==1) {
1003 eANNPair=eANNPairCaseB;
1004 Log(2,"EdbShowerAlg_GS::ReloadANNs","Set eANNPair to eANNPairCaseB");
1005 cout << "EdbShowerAlg_GS::ReloadANNs: eANNPairCaseB at " << eANNPairCaseB << endl;
1006 }
1007 else if (RecoMode==2) {
1008 eANNPair=eANNPairCaseC;
1009 Log(2,"EdbShowerAlg_GS::ReloadANNs","Set eANNPair to eANNPairCaseC");
1010 cout << "EdbShowerAlg_GS::ReloadANNs: eANNPairCaseC at " << eANNPairCaseC << endl;
1011 }
1012 else {
1013 cout << "ELSE RETURN" << endl;
1014 }
1015 Log(2,"EdbShowerAlg_GS::ReloadANNs","ReloadANNs()...done.");
1016}

◆ SelectHighestPInMCArray()

TObjArray * EdbShowerAlg_GS::SelectHighestPInMCArray ( TObjArray *  BTArray)

Helper function for MC purposes:
This function inputs a EdbSegP (Initiator basetrack array) and
for all MC Events it looks for the highest P of the BT.
Mainly used for looking the first BT of a photon/electron shower.
Affects only MC Events, on data it should have no effect.
Only internal usage, should not needed by the everyday user.

1175{
1182
1183 Log(2,"EdbShowerAlg_GS::SelectHighestPInMCArray","SelectHighestPInMCArray()");
1184
1185 if (BTArray->GetEntries()==0) return NULL;
1186
1187 TObjArray* HighestPBTArray = new TObjArray();
1188
1189 Float_t maxP[999999];
1190 Int_t index[999999];
1191 for (int i=0; i<100000; ++i) {
1192 maxP[i]=0;
1193 index[i]=-1;
1194 }
1195
1196 for (int i=0; i<BTArray->GetEntries(); ++i) {
1197 EdbSegP* s = (EdbSegP*)BTArray->At(i);
1198 if (s->MCEvt()<0) continue;
1199 if (s->P() > maxP[s->MCEvt()]) {
1200 maxP[s->MCEvt()] = s->P();
1201 index[s->MCEvt()]=i;
1202 }
1203 }
1204
1205 for (int i=0; i<BTArray->GetEntries(); ++i) {
1206 EdbSegP* s = (EdbSegP*)BTArray->At(i);
1207 if (s->MCEvt()<0) HighestPBTArray->Add(s);
1208 if (s->MCEvt()>0) {
1209 if (s->P()< maxP[s->MCEvt()] ) {
1210 // do nothing, since it is not maximum for this event.
1211 }
1212 else {
1213 HighestPBTArray->Add(s);
1214 }
1215 }
1216 }
1217
1218 if (gEDBDEBUGLEVEL>2) {
1219 cout << "EdbShowerAlg_GS::SelectHighestPInMCArray BTArray->GetEntries() " << BTArray->GetEntries() << endl;
1220 cout << "EdbShowerAlg_GS::SelectHighestPInMCArray HighestPBTArray->GetEntries() " << HighestPBTArray->GetEntries() << endl;
1221
1222 cout << "EdbShowerAlg_GS::SelectHighestPInMCArray Print BTArray: " << BTArray->GetEntries() << endl;
1223
1224 for (int i=0; i<BTArray->GetEntries(); ++i) {
1225 EdbSegP* s = (EdbSegP*)BTArray->At(i);
1226 s->PrintNice();
1227 }
1228 cout << "EdbShowerAlg_GS::SelectHighestPInMCArray Print HighestPBTArray: " << HighestPBTArray->GetEntries() << endl;
1229 for (int i=0; i<HighestPBTArray->GetEntries(); ++i) {
1230 EdbSegP* s = (EdbSegP*)HighestPBTArray->At(i);
1231 s->PrintNice();
1232 }
1233 }
1234
1235 Log(2,"EdbShowerAlg_GS::SelectHighestPInMCArray","SelectHighestPInMCArray()...done.");
1236
1237 return HighestPBTArray;
1238}

◆ Set0()

void EdbShowerAlg_GS::Set0 ( )

Reset all internal variables.

757{
759
760 Log(2,"EdbShowerAlg_GS::EdbShowerAlg_GS","Set0()");
761
762 eRecoMode=0;
763 eCutModeFull=kTRUE;
764
765 eInVtxArrayN=0;
766 eInVtxArray->Clear();
767
768 eInBTArrayN=0;
769 eInBTArray->Clear();
770
771 eRecoShowerArrayN=0;
772 eRecoShowerArray->Clear();
773
774 eAlgName="GS";
775 eAlgValue=999;
776 eANNPairCut[0]=0.5;
777 eANNPairCut[1]=0.48;
778 eANNPairCut[2]=0.45;
779 eSetCleanPairs=kTRUE;
780 eFindPairsPreselected=kTRUE;
781
782 eParaValue[0]=393;
783 eParaString[0]="PARA_GS_CUT_dIP";
784 eParaValue[1]=35;
785 eParaString[1]="PARA_GS_CUT_dMin";
786 eParaValue[2]=85;
787 eParaString[2]="PARA_GS_CUT_dR";
788 eParaValue[3]=25000;
789 eParaString[3]="PARA_GS_CUT_dZ";
790 eParaValue[4]=0.11;
791 eParaString[4]="PARA_GS_CUT_dtheta";
792 eParaValue[5]=1;
793 eParaString[5]="PARA_GS_CUT_PIDDIFF";
794 eParaValue[6]=1;
795 eParaString[6]="PARA_GS_CUT_OPPOSITEFLAG";
796
797
798 eSegmentPIDArray->Reset();
799 eSegmentIDArray->Reset();
800 eSegment2PIDArray->Reset();
801 eSegment2IDArray->Reset();
802
803 Log(2,"EdbShowerAlg_GS::EdbShowerAlg_GS","Set0()...done.");
804 return;
805}

◆ SetCleanPairs()

void EdbShowerAlg_GS::SetCleanPairs ( Bool_t  CleanPairs)
inline
267 {
268 eSetCleanPairs = CleanPairs;
269 }

◆ SetCutModeFull()

void EdbShowerAlg_GS::SetCutModeFull ( Bool_t  CutModeFull)
inline
271 {
272 eCutModeFull = CutModeFull;
273 }

◆ SetFindPairsPreselected()

void EdbShowerAlg_GS::SetFindPairsPreselected ( Bool_t  FindPairsPreselected)
inline
263 {
264 eFindPairsPreselected = FindPairsPreselected;
265 }

◆ SetInVtx()

void EdbShowerAlg_GS::SetInVtx ( EdbVertex *  vtx)
1032{
1033 Log(3,"EdbShowerAlg_GS::SetInVtx","SetInVtx()");
1034 if (eInVtxArrayN!=0) {
1035 Log(2,"EdbShowerAlg_GS::SetInVtx","SetInVtx() WARNING! Array not empty. Clear/Reset it!");
1036 eInVtxArray -> Clear();
1037 eInVtxArrayN=0;
1038 }
1039 eInVtxArray->Add(vtx);
1040 ++eInVtxArrayN;
1041 eInVtxArraySet=kTRUE;
1042 cout << "EdbShowerAlg_GS::SetInVtx Added One Vertex. Now there are " << eInVtxArrayN << " InVtx stored." << endl;
1043 Log(3,"EdbShowerAlg_GS::SetInVtx","SetInVtx()...done");
1044 return;
1045}

◆ SetInVtxArray()

void EdbShowerAlg_GS::SetInVtxArray ( TObjArray *  InVtxArray)
inline
248 {
249 eInVtxArray = InVtxArray;
250 eInVtxArrayN = eInVtxArray->GetEntries();
251 cout << "SetInVtxArray:: " << eInVtxArrayN << " entries set"<<endl;
252 }

◆ SetRecoMode()

void EdbShowerAlg_GS::SetRecoMode ( Int_t  RecoMode)

Manually Set Reco Mode A or B or C
Check if InVtxArray is Filled.
If not, do nothing, cause it does make no sense when
we wanna do Reco Mode and have no vertices filled.

1326{
1331
1332 Log(2,"EdbShowerAlg_GS::SetRecoMode","SetRecoMode()...");
1333
1334 if (eRecoMode>2) {
1335 Log(2,"EdbShowerAlg_GS::SetRecoMode","WARNING eRecoMode>2!!!");
1336 Log(2,"EdbShowerAlg_GS::SetRecoMode","Set automatically Reco Mode C (eRecoMode=2).");
1337 eRecoMode=2;
1338 }
1339
1340 if (eInVtxArrayN==0) {
1341
1342 // So, welches nehmen wir nun denn jetzt by default??
1343 // Die ROC Curve sagt eigentlich, dass CASE C besser sei...
1344
1345 Log(2,"EdbShowerAlg_GS::SetRecoMode","No InVtxArray there. Will not set Reco Mode A.");
1346 Log(2,"EdbShowerAlg_GS::SetRecoMode","Set automatically Reco Mode B (eRecoMode=1).");
1347 eRecoMode=1;
1348
1349 Log(2,"EdbShowerAlg_GS::SetRecoMode","No InVtxArray there. Will not set Reco Mode A.");
1350 Log(2,"EdbShowerAlg_GS::SetRecoMode","Set automatically Reco Mode C (eRecoMode=2).");
1351 eRecoMode=2;
1352
1353 }
1354 else {
1355 eRecoMode=RecoMode;
1356 ReloadANNs(RecoMode);
1357 }
1358
1359 cout << "DEBUG eRecoMode " << eRecoMode << endl;
1360
1361 // Set the dZ cutvalue to just one plate after, since we start
1362 // in these cases from real Initiator Basetracks (and not vertices)
1363 // so we just look for pairs directly related to this InBT.
1364 if (eRecoMode==1 || eRecoMode == 2) {
1365 SetParameter(3,1800);
1366 }
1367 else {
1368 SetParameter(3,25000);
1369 }
1370
1371 Log(2,"EdbShowerAlg_GS::SetRecoMode","SetRecoMode() eRecoMode = %d ",eRecoMode);
1372 Log(2,"EdbShowerAlg_GS::SetRecoMode","SetRecoMode()...done.");
1373 return;
1374}
EdbShowerAlg::SetParameter
void SetParameter(Int_t parNr, Float_t parvalue)
Definition: EdbShowerAlg.cxx:82

Member Data Documentation

◆ eANNPair

TMultiLayerPerceptron* EdbShowerAlg_GS::eANNPair
private

◆ eANNPairCaseA

TMultiLayerPerceptron* EdbShowerAlg_GS::eANNPairCaseA
private

◆ eANNPairCaseB

TMultiLayerPerceptron* EdbShowerAlg_GS::eANNPairCaseB
private

◆ eANNPairCaseC

TMultiLayerPerceptron* EdbShowerAlg_GS::eANNPairCaseC
private

◆ eANNPairCut

Float_t EdbShowerAlg_GS::eANNPairCut[3]
private

◆ eANNPairTree

TTree* EdbShowerAlg_GS::eANNPairTree
private

◆ eCutModeFull

Bool_t EdbShowerAlg_GS::eCutModeFull
private

◆ eFindPairsPreselected

Bool_t EdbShowerAlg_GS::eFindPairsPreselected
private

◆ eInVtxArray

TObjArray* EdbShowerAlg_GS::eInVtxArray
private

◆ eInVtxArrayN

Int_t EdbShowerAlg_GS::eInVtxArrayN
private

◆ eInVtxArraySet

Bool_t EdbShowerAlg_GS::eInVtxArraySet
private

◆ eRecoMode

Int_t EdbShowerAlg_GS::eRecoMode
private

◆ eSegment2IDArray

TArrayI* EdbShowerAlg_GS::eSegment2IDArray
private

◆ eSegment2PIDArray

TArrayI* EdbShowerAlg_GS::eSegment2PIDArray
private

◆ eSegmentIDArray

TArrayI* EdbShowerAlg_GS::eSegmentIDArray
private

◆ eSegmentPIDArray

TArrayI* EdbShowerAlg_GS::eSegmentPIDArray
private

◆ eSetCleanPairs

Bool_t EdbShowerAlg_GS::eSetCleanPairs
private

◆ eValueGSNN_var00

Float_t EdbShowerAlg_GS::eValueGSNN_var00
private

◆ eValueGSNN_var01

Float_t EdbShowerAlg_GS::eValueGSNN_var01
private

◆ eValueGSNN_var02

Float_t EdbShowerAlg_GS::eValueGSNN_var02
private

◆ eValueGSNN_var03

Float_t EdbShowerAlg_GS::eValueGSNN_var03
private

◆ eValueGSNN_var04

Float_t EdbShowerAlg_GS::eValueGSNN_var04
private

◆ eValueGSNN_var05

Float_t EdbShowerAlg_GS::eValueGSNN_var05
private

◆ eValueGSNN_var06

Float_t EdbShowerAlg_GS::eValueGSNN_var06
private

◆ eValueGSNN_varInput

Float_t EdbShowerAlg_GS::eValueGSNN_varInput
private

◆ eValueGSNN_varOutput

Float_t EdbShowerAlg_GS::eValueGSNN_varOutput
private
The documentation for this class was generated from the following files: