From c502f100ebf5d4eb35b465f5054fe6f1c685bbf2 Mon Sep 17 00:00:00 2001 From: Fabio Date: Thu, 27 Feb 2025 18:36:43 +0100 Subject: [PATCH 01/19] create new tree with flag0 only --- src/appl/emrec/emvertex.cpp | 17 +++++++++++++---- 1 file changed, 13 insertions(+), 4 deletions(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index 9a74be75..17927d64 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -229,6 +229,8 @@ void ReadVertex(EdbID id, TEnv &env) gEVR.eQualityMode= env.GetValue("emvertex.vtx.QualityMode" , 0); // (0:=Prob/(sigVX^2+sigVY^2); 1:= inverse average track-vertex distance) TCut cutvtx = env.GetValue("emvertex.vtx.cutvtx" , "(flag==0||flag==3)&&n>4"); + TObjArray* v_out = new TObjArray(); + EdbDataProc *dproc = new EdbDataProc(); TString name; gSproc.MakeFileName(name,id,"vtx.root",false); @@ -241,7 +243,8 @@ void ReadVertex(EdbID id, TEnv &env) EdbPVRec *vtr = new EdbPVRec(); vtr->SetScanCond( new EdbScanCond(gCond) ); gSproc.ReadTracksTree( idset,*vtr, cuttr); - AddCompatibleTracks( *vtr, gAli ); // assign to the vertices of gAli additional tracks from vtr if any + AddCompatibleTracks( *vtr, gAli , v_out); // assign to the vertices of gAli additional tracks from vtr if any + EdbDataProc::MakeVertexTree(v_out,"flag0.vtx.root"); } } } @@ -315,21 +318,27 @@ void MakeScanCondBT(EdbScanCond &cond, TEnv &env) cond.SetName("SND_basetrack"); } -void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx) +void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out) { int ntr = v_trk.Ntracks(); int nvtx = v_vtx.Nvtx(); Log(1,"AddCompatibleTracks", "%d tracks, %d vertex", ntr,nvtx ); for(int iv=0; ivSetFlag(999999); - v_vtx.AddTrack(t); + flag1 = true; + t->SetFlag(999999); + v_vtx.AddTrack(t); + } } + if (!flag1) { + v_out.AddVertex(v); + } } } } From 0595167d92031d811c0d3e8cf421b7c4f9836167 Mon Sep 17 00:00:00 2001 From: Fabio Date: Thu, 27 Feb 2025 18:43:02 +0100 Subject: [PATCH 02/19] fix tobj --- src/appl/emrec/emvertex.cpp | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index 17927d64..0f967030 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -244,7 +244,7 @@ void ReadVertex(EdbID id, TEnv &env) vtr->SetScanCond( new EdbScanCond(gCond) ); gSproc.ReadTracksTree( idset,*vtr, cuttr); AddCompatibleTracks( *vtr, gAli , v_out); // assign to the vertices of gAli additional tracks from vtr if any - EdbDataProc::MakeVertexTree(v_out,"flag0.vtx.root"); + EdbDataProc::MakeVertexTree(&v_out,"flag0.vtx.root"); } } } @@ -337,8 +337,7 @@ void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out) } } if (!flag1) { - v_out.AddVertex(v); - } + v_out.Add(v); } } } From 1f04915efa803905e4e15cf5e60364e333e1e473 Mon Sep 17 00:00:00 2001 From: Fabio Date: Thu, 27 Feb 2025 18:44:28 +0100 Subject: [PATCH 03/19] fix declaration --- src/appl/emrec/emvertex.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index 0f967030..9c815318 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -26,7 +26,7 @@ void ReadVertex(EdbID id,TEnv &env); void MakeScanCondBT(EdbScanCond &cond, TEnv &env); void SetTracksErrors(TObjArray &tracks, EdbScanCond &cond, float p, float m); void do_vertex(TEnv &env); -void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx); +void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out); bool IsCompatible(EdbVertex &v, EdbTrackP &t); void Display( const char *dsname, EdbVertexRec *evr, TEnv &env ); From 102314ef3a59dd2e92fc699c188a98fdb244937e Mon Sep 17 00:00:00 2001 From: Fabio Date: Thu, 27 Feb 2025 19:11:40 +0100 Subject: [PATCH 04/19] exclude vtx tracks --- src/appl/emrec/emvertex.cpp | 12 ++++++++++-- 1 file changed, 10 insertions(+), 2 deletions(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index 9c815318..4b928dd4 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -229,7 +229,7 @@ void ReadVertex(EdbID id, TEnv &env) gEVR.eQualityMode= env.GetValue("emvertex.vtx.QualityMode" , 0); // (0:=Prob/(sigVX^2+sigVY^2); 1:= inverse average track-vertex distance) TCut cutvtx = env.GetValue("emvertex.vtx.cutvtx" , "(flag==0||flag==3)&&n>4"); - TObjArray* v_out = new TObjArray(); + TObjArray v_out; EdbDataProc *dproc = new EdbDataProc(); TString name; @@ -244,7 +244,7 @@ void ReadVertex(EdbID id, TEnv &env) vtr->SetScanCond( new EdbScanCond(gCond) ); gSproc.ReadTracksTree( idset,*vtr, cuttr); AddCompatibleTracks( *vtr, gAli , v_out); // assign to the vertices of gAli additional tracks from vtr if any - EdbDataProc::MakeVertexTree(&v_out,"flag0.vtx.root"); + EdbDataProc::MakeVertexTree(v_out,"flag0.vtx.root"); } } } @@ -326,10 +326,18 @@ void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out) for(int iv=0; iv trackids; EdbVertex *v = v_vtx.GetVertex(iv); + for(int i=0; iN(); i++){ + EdbTrackP *t = (EdbTrackP*)v->GetTrack(i); + int trid = t->ID(); + trackids.push_back(trid); + } for(int it=0; itID(); + if (std::find(trackids.begin(), trackids.end(), trid)!=trackids.end()) continue; if( IsCompatible(*v,*t) ) { flag1 = true; t->SetFlag(999999); From 767f60feb23928c1ab4d6b2a1ff8007e9c0c0742 Mon Sep 17 00:00:00 2001 From: Fabio Date: Mon, 17 Mar 2025 16:40:12 +0100 Subject: [PATCH 05/19] take only upstream tracks, only closest track --- src/appl/emrec/emvertex.cpp | 30 ++++++++++++++++++++---------- src/libEdr/EdbPattern.cxx | 12 +++++++----- 2 files changed, 27 insertions(+), 15 deletions(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index 4b928dd4..e72a1e27 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -26,7 +26,7 @@ void ReadVertex(EdbID id,TEnv &env); void MakeScanCondBT(EdbScanCond &cond, TEnv &env); void SetTracksErrors(TObjArray &tracks, EdbScanCond &cond, float p, float m); void do_vertex(TEnv &env); -void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out); +void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TH1F* h_r2, TH1F* h_dz); bool IsCompatible(EdbVertex &v, EdbTrackP &t); void Display( const char *dsname, EdbVertexRec *evr, TEnv &env ); @@ -243,7 +243,9 @@ void ReadVertex(EdbID id, TEnv &env) EdbPVRec *vtr = new EdbPVRec(); vtr->SetScanCond( new EdbScanCond(gCond) ); gSproc.ReadTracksTree( idset,*vtr, cuttr); - AddCompatibleTracks( *vtr, gAli , v_out); // assign to the vertices of gAli additional tracks from vtr if any + TH1F *h_r2 = new TH1F("h_r2", "h_r2", 50,0,5); + TH1F *h_dz = new TH1F("h_dz", "h_dz", 400,0,4000); + AddCompatibleTracks( *vtr, gAli , v_out, h_r2, h_dz); // assign to the vertices of gAli additional tracks from vtr if any EdbDataProc::MakeVertexTree(v_out,"flag0.vtx.root"); } } @@ -318,7 +320,7 @@ void MakeScanCondBT(EdbScanCond &cond, TEnv &env) cond.SetName("SND_basetrack"); } -void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out) +void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TH1F* h_r2, TH1F* h_dz) { int ntr = v_trk.Ntracks(); int nvtx = v_vtx.Nvtx(); @@ -333,32 +335,40 @@ void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out) int trid = t->ID(); trackids.push_back(trid); } + EdbTrackP *t_chosen = 0; + float dr2 = 1e9f; for(int it=0; itID(); if (std::find(trackids.begin(), trackids.end(), trid)!=trackids.end()) continue; - if( IsCompatible(*v,*t) ) { + if( IsCompatible(*v,*t, &r2, &dz) ) { flag1 = true; t->SetFlag(999999); - v_vtx.AddTrack(t); - } + if( r2 < dr2 ) { t_chosen=t; } } + } + h_r2->Fill(r2); + h_dz->Fill(dz); + v_vtx.AddTrack(t); if (!flag1) { v_out.Add(v); } } } -bool IsCompatible(EdbVertex &v, EdbTrackP &t) +bool IsCompatible(EdbVertex &v, EdbTrackP &t, float *r2, float *dz) { EdbSegP ss; + EdbSegP tst = t.GetSegmentFirst(); + float tz = tst.Z(); + if (tz > v.VZ()) return false; //only tracks starting upstream of the vtx t.EstimatePositionAt(v.VZ(),ss); float dx=ss.X()-v.VX(); float dy=ss.Y()-v.VY(); - float r2 = Sqrt(dx*dx+dy*dy); - float dz = Abs(ss.DZ()); - if(r2<5&&dz<4000) { printf("r2=%.4f dz=%.2f\n",r2,ss.DZ()); return true;} + *r2 = Sqrt(dx*dx+dy*dy); + *dz = Abs(ss.DZ()); + if(*r2<5&&*dz<4000) { printf("r2=%.4f dz=%.2f\n",*r2,ss.DZ()); return true;} return false; } diff --git a/src/libEdr/EdbPattern.cxx b/src/libEdr/EdbPattern.cxx index 6218ab33..ef10f07f 100644 --- a/src/libEdr/EdbPattern.cxx +++ b/src/libEdr/EdbPattern.cxx @@ -1096,7 +1096,7 @@ int EdbTrackP::EstimatePositionAt( Float_t z, EdbSegP &ss ) // TODO: dz=0: mean? if( N()<2 ) return 0; EdbSegP *s1=0,*s2=0; - float dz1,dz2; dz1=dz2=kMaxInt; + float dz1,dz2; dz1=dz2=1e9f; //kMaxInt into float; for(int i=0; iX()-s1->X(); float dy0 = s2->Y()-s1->Y(); float dz = z-s1->Z(); - ss.SetX( s1->X() + dz*dx0/dz0 ); - ss.SetY( s1->Y() + dz*dy0/dz0 ); + float tx = dx0/dz0; + float ty = dy0/dz0; + ss.SetX( s1->X() + dz*tx ); + ss.SetY( s1->Y() + dz*ty ); ss.SetZ( z ); ss.SetDZ( dz ); // keep dz distance - ss.SetTX( dx0/dz0 ); - ss.SetTY( dy0/dz0 ); + ss.SetTX( tx ); + ss.SetTY( ty ); return 1; } From 45f898e8a589e0d47ba02a2d486fbf8df8ca3ebd Mon Sep 17 00:00:00 2001 From: Fabio Date: Mon, 17 Mar 2025 16:57:21 +0100 Subject: [PATCH 06/19] clean logs --- src/appl/emrec/emvertex.cpp | 17 +++++++++-------- src/libEdr/EdbVertex.cxx | 2 +- 2 files changed, 10 insertions(+), 9 deletions(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index e72a1e27..2105f531 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -27,7 +27,7 @@ void MakeScanCondBT(EdbScanCond &cond, TEnv &env); void SetTracksErrors(TObjArray &tracks, EdbScanCond &cond, float p, float m); void do_vertex(TEnv &env); void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TH1F* h_r2, TH1F* h_dz); -bool IsCompatible(EdbVertex &v, EdbTrackP &t); +bool IsCompatible(EdbVertex &v, EdbTrackP &t, float *r2, float *dz); void Display( const char *dsname, EdbVertexRec *evr, TEnv &env ); //---------------------------------------------------------------------------------------- @@ -336,21 +336,22 @@ void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TH1 trackids.push_back(trid); } EdbTrackP *t_chosen = 0; - float dr2 = 1e9f; + float r2, dz, r2max; r2max=1e9f; for(int it=0; itID(); if (std::find(trackids.begin(), trackids.end(), trid)!=trackids.end()) continue; if( IsCompatible(*v,*t, &r2, &dz) ) { - flag1 = true; - t->SetFlag(999999); - if( r2 < dr2 ) { t_chosen=t; } + flag1 = true; + t->SetFlag(999999); + if( r2 < r2max ) { t_chosen=t; } } } h_r2->Fill(r2); h_dz->Fill(dz); - v_vtx.AddTrack(t); + v_vtx.AddTrack(t_chosen); + Log(1,"AddCompatibleTracks","Closest track found at r2=%.4f dz=%.2f\n",r2,dz); if (!flag1) { v_out.Add(v); } @@ -360,8 +361,8 @@ void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TH1 bool IsCompatible(EdbVertex &v, EdbTrackP &t, float *r2, float *dz) { EdbSegP ss; - EdbSegP tst = t.GetSegmentFirst(); - float tz = tst.Z(); + EdbSegP *tst = t.GetSegmentFirst(); + float tz = tst->Z(); if (tz > v.VZ()) return false; //only tracks starting upstream of the vtx t.EstimatePositionAt(v.VZ(),ss); float dx=ss.X()-v.VX(); diff --git a/src/libEdr/EdbVertex.cxx b/src/libEdr/EdbVertex.cxx index 4813c913..b8dd793d 100644 --- a/src/libEdr/EdbVertex.cxx +++ b/src/libEdr/EdbVertex.cxx @@ -160,7 +160,7 @@ EdbSegP *EdbVertex::GetTrackV(int i, bool usesegpar) { EdbTrackP *t = GetTrack(i); if(!t) return 0; EdbSegP *s = t->TrackExtremity(Zpos(i), usesegpar); - if( s->P()>=0 && (s->P() != t->P()) ) Log(1,"GetTrackV","Warning! segment momentum=%f is not equal to the track momentum=%f",s->P(), t->P() ); + //if( s->P()>=0 && (s->P() != t->P()) ) Log(1,"GetTrackV","Warning! segment momentum=%f is not equal to the track momentum=%f",s->P(), t->P() ); return s; } From 36ae09ae6405a9d4f61208caf71f31062a19c532 Mon Sep 17 00:00:00 2001 From: Fabio Date: Mon, 17 Mar 2025 17:53:31 +0100 Subject: [PATCH 07/19] add par to env --- src/appl/emrec/emvertex.cpp | 37 +++++++++++++++++++++++-------------- 1 file changed, 23 insertions(+), 14 deletions(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index 2105f531..e7dbd206 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -26,8 +26,8 @@ void ReadVertex(EdbID id,TEnv &env); void MakeScanCondBT(EdbScanCond &cond, TEnv &env); void SetTracksErrors(TObjArray &tracks, EdbScanCond &cond, float p, float m); void do_vertex(TEnv &env); -void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TH1F* h_r2, TH1F* h_dz); -bool IsCompatible(EdbVertex &v, EdbTrackP &t, float *r2, float *dz); +void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TH1F* h_r2, TH1F* h_dz); +bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz); void Display( const char *dsname, EdbVertexRec *evr, TEnv &env ); //---------------------------------------------------------------------------------------- @@ -68,6 +68,9 @@ void set_default(TEnv &env) env.SetValue("emvertex.trfit.doit" , 1 ); env.SetValue("emvertex.trfit.P" , 10 ); env.SetValue("emvertex.trfit.M" , 0.139); + env.SetValue("emvertex.trfit.r2max", 5. ); + env.SetValue("emvertex.trfit.dzmax", 4000. ); + env.SetValue("emvertex.bt.Sigma0", "0.2 0.2 0.002 0.002" ); env.SetValue("emvertex.bt.Degrad", 5. ); } @@ -245,7 +248,7 @@ void ReadVertex(EdbID id, TEnv &env) gSproc.ReadTracksTree( idset,*vtr, cuttr); TH1F *h_r2 = new TH1F("h_r2", "h_r2", 50,0,5); TH1F *h_dz = new TH1F("h_dz", "h_dz", 400,0,4000); - AddCompatibleTracks( *vtr, gAli , v_out, h_r2, h_dz); // assign to the vertices of gAli additional tracks from vtr if any + AddCompatibleTracks(env, *vtr, gAli , v_out, h_r2, h_dz); // assign to the vertices of gAli additional tracks from vtr if any EdbDataProc::MakeVertexTree(v_out,"flag0.vtx.root"); } } @@ -320,7 +323,7 @@ void MakeScanCondBT(EdbScanCond &cond, TEnv &env) cond.SetName("SND_basetrack"); } -void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TH1F* h_r2, TH1F* h_dz) +void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TH1F* h_r2, TH1F* h_dz) { int ntr = v_trk.Ntracks(); int nvtx = v_vtx.Nvtx(); @@ -330,35 +333,39 @@ void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TH1 bool flag1 = false; std::vector trackids; EdbVertex *v = v_vtx.GetVertex(iv); + Log(1,"AddCompatibleTracks","Looking for parent tracks of vtx: %i\n",v->ID()); for(int i=0; iN(); i++){ EdbTrackP *t = (EdbTrackP*)v->GetTrack(i); int trid = t->ID(); trackids.push_back(trid); } EdbTrackP *t_chosen = 0; - float r2, dz, r2max; r2max=1e9f; + float r2, dz, r2max, dzmax; r2max=dzmax=1e9f; for(int it=0; itID(); if (std::find(trackids.begin(), trackids.end(), trid)!=trackids.end()) continue; - if( IsCompatible(*v,*t, &r2, &dz) ) { + if( IsCompatible(env, *v,*t, &r2, &dz) ) { flag1 = true; t->SetFlag(999999); - if( r2 < r2max ) { t_chosen=t; } + if( r2 < r2max ) { r2max=r2; dzmax=dz; t_chosen=t; } + // v_vtx.AddTrack(t); } } - h_r2->Fill(r2); - h_dz->Fill(dz); - v_vtx.AddTrack(t_chosen); - Log(1,"AddCompatibleTracks","Closest track found at r2=%.4f dz=%.2f\n",r2,dz); - if (!flag1) { + if (flag1){ + h_r2->Fill(r2); + h_dz->Fill(dz); + v_vtx.AddTrack(t_chosen); + Log(1,"AddCompatibleTracks","Closest track found at r2=%.4f dz=%.2f\n",r2max,dzmax); + } + else { v_out.Add(v); } } } -bool IsCompatible(EdbVertex &v, EdbTrackP &t, float *r2, float *dz) +bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz) { EdbSegP ss; EdbSegP *tst = t.GetSegmentFirst(); @@ -369,7 +376,9 @@ bool IsCompatible(EdbVertex &v, EdbTrackP &t, float *r2, float *dz) float dy=ss.Y()-v.VY(); *r2 = Sqrt(dx*dx+dy*dy); *dz = Abs(ss.DZ()); - if(*r2<5&&*dz<4000) { printf("r2=%.4f dz=%.2f\n",*r2,ss.DZ()); return true;} + float r2max = env.GetValue("emvertex.trfit.r2max" , 5. ); + float dzmax = env.GetValue("emvertex.trfit.dzmax" , 4000. ); + if(*r2 Date: Tue, 18 Mar 2025 17:40:03 +0100 Subject: [PATCH 08/19] add output file with tracks info --- src/appl/emrec/emvertex.cpp | 19 ++++++++++++------- 1 file changed, 12 insertions(+), 7 deletions(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index e7dbd206..aa2b0b91 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -26,7 +26,7 @@ void ReadVertex(EdbID id,TEnv &env); void MakeScanCondBT(EdbScanCond &cond, TEnv &env); void SetTracksErrors(TObjArray &tracks, EdbScanCond &cond, float p, float m); void do_vertex(TEnv &env); -void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TH1F* h_r2, TH1F* h_dz); +void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TNtuple* outTree); bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz); void Display( const char *dsname, EdbVertexRec *evr, TEnv &env ); @@ -246,10 +246,13 @@ void ReadVertex(EdbID id, TEnv &env) EdbPVRec *vtr = new EdbPVRec(); vtr->SetScanCond( new EdbScanCond(gCond) ); gSproc.ReadTracksTree( idset,*vtr, cuttr); - TH1F *h_r2 = new TH1F("h_r2", "h_r2", 50,0,5); - TH1F *h_dz = new TH1F("h_dz", "h_dz", 400,0,4000); - AddCompatibleTracks(env, *vtr, gAli , v_out, h_r2, h_dz); // assign to the vertices of gAli additional tracks from vtr if any + TNtuple *outTree = new TNtuple("tracks","Tree of matched tracks","chosen:n:vid:tid:nseg:npl:tx:ty:firstp:lastp:r2:dz"); + AddCompatibleTracks(env, *vtr, gAli , v_out, outTree); // assign to the vertices of gAli additional tracks from vtr if any EdbDataProc::MakeVertexTree(v_out,"flag0.vtx.root"); + TFile *outFile = new TFile("found_tracks.root","RECREATE"); + outTree->Write(); + outFile->Write(); + outFile->Close(); } } } @@ -323,7 +326,7 @@ void MakeScanCondBT(EdbScanCond &cond, TEnv &env) cond.SetName("SND_basetrack"); } -void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TH1F* h_r2, TH1F* h_dz) +void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TNtuple* outTree) { int ntr = v_trk.Ntracks(); int nvtx = v_vtx.Nvtx(); @@ -341,6 +344,7 @@ void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray } EdbTrackP *t_chosen = 0; float r2, dz, r2max, dzmax; r2max=dzmax=1e9f; + int founds=0; for(int it=0; itSetFlag(999999); if( r2 < r2max ) { r2max=r2; dzmax=dz; t_chosen=t; } // v_vtx.AddTrack(t); + founds++; + outTree->Fill(0, 1, v->ID(), t->ID() ,t->N(), t->Npl(), t->TX(), t->TY(), t->GetSegmentFirst()->Plate(), t->GetSegmentLast()->Plate(), r2, dz); } } if (flag1){ - h_r2->Fill(r2); - h_dz->Fill(dz); v_vtx.AddTrack(t_chosen); + outTree->Fill(1, founds, v->ID(), t_chosen->ID(), t_chosen->N(), t_chosen->Npl(), t_chosen->TX(), t_chosen->TY(), t_chosen->GetSegmentFirst()->Plate(), t_chosen->GetSegmentLast()->Plate(), r2max, dzmax); Log(1,"AddCompatibleTracks","Closest track found at r2=%.4f dz=%.2f\n",r2max,dzmax); } else { From e783e1c211ab479747e94f326977f0d1e06abbbd Mon Sep 17 00:00:00 2001 From: dcentanni Date: Wed, 19 Mar 2025 10:51:16 +0100 Subject: [PATCH 09/19] Vertex files with flag0 and flag1 created + tunable trackfind --- src/appl/emrec/emvertex.cpp | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index aa2b0b91..a2bc65df 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -26,7 +26,7 @@ void ReadVertex(EdbID id,TEnv &env); void MakeScanCondBT(EdbScanCond &cond, TEnv &env); void SetTracksErrors(TObjArray &tracks, EdbScanCond &cond, float p, float m); void do_vertex(TEnv &env); -void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TNtuple* outTree); +void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TObjArray &v_out2, TNtuple* outTree); bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz); void Display( const char *dsname, EdbVertexRec *evr, TEnv &env ); @@ -233,6 +233,7 @@ void ReadVertex(EdbID id, TEnv &env) TCut cutvtx = env.GetValue("emvertex.vtx.cutvtx" , "(flag==0||flag==3)&&n>4"); TObjArray v_out; + TObjArray v_out2; EdbDataProc *dproc = new EdbDataProc(); TString name; @@ -247,8 +248,9 @@ void ReadVertex(EdbID id, TEnv &env) vtr->SetScanCond( new EdbScanCond(gCond) ); gSproc.ReadTracksTree( idset,*vtr, cuttr); TNtuple *outTree = new TNtuple("tracks","Tree of matched tracks","chosen:n:vid:tid:nseg:npl:tx:ty:firstp:lastp:r2:dz"); - AddCompatibleTracks(env, *vtr, gAli , v_out, outTree); // assign to the vertices of gAli additional tracks from vtr if any + AddCompatibleTracks(env, *vtr, gAli , v_out, , v_out2, outTree); // assign to the vertices of gAli additional tracks from vtr if any EdbDataProc::MakeVertexTree(v_out,"flag0.vtx.root"); + EdbDataProc::MakeVertexTree(v_out2,"flag1.vtx.root"); TFile *outFile = new TFile("found_tracks.root","RECREATE"); outTree->Write(); outFile->Write(); @@ -326,7 +328,7 @@ void MakeScanCondBT(EdbScanCond &cond, TEnv &env) cond.SetName("SND_basetrack"); } -void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TNtuple* outTree) +void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TObjArray &v_out2, TNtuple* outTree) { int ntr = v_trk.Ntracks(); int nvtx = v_vtx.Nvtx(); @@ -367,6 +369,7 @@ void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray else { v_out.Add(v); } + else { v_out2.Add(v);} } } From 0d3ce43cbad53f8d99292247aecbe37a10a7bb76 Mon Sep 17 00:00:00 2001 From: dcentanni Date: Mon, 24 Mar 2025 14:48:38 +0100 Subject: [PATCH 10/19] Some fixes --- src/appl/emrec/emvertex.cpp | 11 ++++++++--- 1 file changed, 8 insertions(+), 3 deletions(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index a2bc65df..5a4ecc01 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -167,6 +167,7 @@ int main(int argc, char* argv[]) bool do_set = false; bool do_display = false; bool do_read = false; + bool do_vtxrefit = false; for(int i=1; iSetScanCond( new EdbScanCond(gCond) ); gSproc.ReadTracksTree( idset,*vtr, cuttr); TNtuple *outTree = new TNtuple("tracks","Tree of matched tracks","chosen:n:vid:tid:nseg:npl:tx:ty:firstp:lastp:r2:dz"); - AddCompatibleTracks(env, *vtr, gAli , v_out, , v_out2, outTree); // assign to the vertices of gAli additional tracks from vtr if any + AddCompatibleTracks(env, *vtr, gAli , v_out, v_out2, outTree); // assign to the vertices of gAli additional tracks from vtr if any EdbDataProc::MakeVertexTree(v_out,"flag0.vtx.root"); EdbDataProc::MakeVertexTree(v_out2,"flag1.vtx.root"); TFile *outFile = new TFile("found_tracks.root","RECREATE"); @@ -351,7 +356,7 @@ void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray { EdbTrackP *t = v_trk.GetTrack(it); int trid = t->ID(); - if (std::find(trackids.begin(), trackids.end(), trid)!=trackids.end()) continue; + if (std::find(trackids.begin(), trackids.end(), trid)!=trackids.end()) continue; //Maybe here can be changed to EdbVertex::TrackInVertex(EdbTrackP *t) if( IsCompatible(env, *v,*t, &r2, &dz) ) { flag1 = true; t->SetFlag(999999); @@ -365,11 +370,11 @@ void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray v_vtx.AddTrack(t_chosen); outTree->Fill(1, founds, v->ID(), t_chosen->ID(), t_chosen->N(), t_chosen->Npl(), t_chosen->TX(), t_chosen->TY(), t_chosen->GetSegmentFirst()->Plate(), t_chosen->GetSegmentLast()->Plate(), r2max, dzmax); Log(1,"AddCompatibleTracks","Closest track found at r2=%.4f dz=%.2f\n",r2max,dzmax); + v_out2.Add(v); } else { v_out.Add(v); } - else { v_out2.Add(v);} } } From 7f7f8c3c6af1fbe2987c2828fad4b4a3ed0658e2 Mon Sep 17 00:00:00 2001 From: Fabio Date: Fri, 21 Mar 2025 16:47:51 +0100 Subject: [PATCH 11/19] add imp --- src/appl/emrec/emvertex.cpp | 32 +++++++++++++++++++++++++------- 1 file changed, 25 insertions(+), 7 deletions(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index 5a4ecc01..702966f7 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -29,6 +29,7 @@ void do_vertex(TEnv &env); void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TObjArray &v_out2, TNtuple* outTree); bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz); void Display( const char *dsname, EdbVertexRec *evr, TEnv &env ); +float distance(const EdbTrackP& t, const EdbVertex& v); //---------------------------------------------------------------------------------------- void print_help_message() @@ -70,6 +71,7 @@ void set_default(TEnv &env) env.SetValue("emvertex.trfit.M" , 0.139); env.SetValue("emvertex.trfit.r2max", 5. ); env.SetValue("emvertex.trfit.dzmax", 4000. ); + env.SetValue("emvertex.trfit.impmax", 5. ); env.SetValue("emvertex.bt.Sigma0", "0.2 0.2 0.002 0.002" ); env.SetValue("emvertex.bt.Degrad", 5. ); @@ -260,6 +262,8 @@ void ReadVertex(EdbID id, TEnv &env) outTree->Write(); outFile->Write(); outFile->Close(); + delete outTree; + delete outFile; } } } @@ -350,7 +354,7 @@ void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray trackids.push_back(trid); } EdbTrackP *t_chosen = 0; - float r2, dz, r2max, dzmax; r2max=dzmax=1e9f; + float r2, dz, imp, imp2, r2max, dzmax, impmax, impmax2; r2max=dzmax=impmax=impmax2=1e9f; int founds=0; for(int it=0; itSetFlag(999999); - if( r2 < r2max ) { r2max=r2; dzmax=dz; t_chosen=t; } + if( r2 < r2max ) { r2max=r2; dzmax=dz; impmax=imp; impmax2=imp2; t_chosen=t; } // v_vtx.AddTrack(t); founds++; - outTree->Fill(0, 1, v->ID(), t->ID() ,t->N(), t->Npl(), t->TX(), t->TY(), t->GetSegmentFirst()->Plate(), t->GetSegmentLast()->Plate(), r2, dz); + outTree->Fill(0, 1, v->ID(), t->ID() ,t->N(), t->Npl(), t->TX(), t->TY(), t->GetSegmentFirst()->Plate(), t->GetSegmentLast()->Plate(), r2, dz, imp, imp2); } } if (flag1){ v_vtx.AddTrack(t_chosen); - outTree->Fill(1, founds, v->ID(), t_chosen->ID(), t_chosen->N(), t_chosen->Npl(), t_chosen->TX(), t_chosen->TY(), t_chosen->GetSegmentFirst()->Plate(), t_chosen->GetSegmentLast()->Plate(), r2max, dzmax); + outTree->Fill(1, founds, v->ID(), t_chosen->ID(), t_chosen->N(), t_chosen->Npl(), t_chosen->TX(), t_chosen->TY(), t_chosen->GetSegmentFirst()->Plate(), t_chosen->GetSegmentLast()->Plate(), r2max, dzmax, impmax, impmax2); Log(1,"AddCompatibleTracks","Closest track found at r2=%.4f dz=%.2f\n",r2max,dzmax); v_out2.Add(v); } @@ -378,7 +382,7 @@ void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray } } -bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz) +bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz, float *imp, float *imp2) { EdbSegP ss; EdbSegP *tst = t.GetSegmentFirst(); @@ -388,10 +392,13 @@ bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz) float dx=ss.X()-v.VX(); float dy=ss.Y()-v.VY(); *r2 = Sqrt(dx*dx+dy*dy); - *dz = Abs(ss.DZ()); + *dz = ss.DZ(); + *imp = distance(t, v); + *imp2 = distance(ss, v); float r2max = env.GetValue("emvertex.trfit.r2max" , 5. ); float dzmax = env.GetValue("emvertex.trfit.dzmax" , 4000. ); - if(*r2FitTrackKFS(); } } + +float distance(const EdbTrackP& t, const EdbVertex& v) { + float dx = v.VX() - t.X(); + float dy = v.VY() - t.Y(); + float dz = v.VZ() - t.Z(); + float tx = t.TX(); + float ty = t.TY(); + float nom = std::pow((dx*ty - dy*tx),2) + std::pow((dy - dz*ty),2) + std::pow((dz*tx - dx),2); + float denom = std::pow((tx),2) + std::pow((ty),2) + 1.; +return sqrt(nom/denom); +} From 1116d8e4441607c5788e728f7a89f886d9e49a87 Mon Sep 17 00:00:00 2001 From: Fabio Date: Fri, 21 Mar 2025 16:50:14 +0100 Subject: [PATCH 12/19] rm imp --- src/appl/emrec/emvertex.cpp | 15 +++++---------- 1 file changed, 5 insertions(+), 10 deletions(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index 702966f7..042211d5 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -29,7 +29,6 @@ void do_vertex(TEnv &env); void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TObjArray &v_out2, TNtuple* outTree); bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz); void Display( const char *dsname, EdbVertexRec *evr, TEnv &env ); -float distance(const EdbTrackP& t, const EdbVertex& v); //---------------------------------------------------------------------------------------- void print_help_message() @@ -71,7 +70,6 @@ void set_default(TEnv &env) env.SetValue("emvertex.trfit.M" , 0.139); env.SetValue("emvertex.trfit.r2max", 5. ); env.SetValue("emvertex.trfit.dzmax", 4000. ); - env.SetValue("emvertex.trfit.impmax", 5. ); env.SetValue("emvertex.bt.Sigma0", "0.2 0.2 0.002 0.002" ); env.SetValue("emvertex.bt.Degrad", 5. ); @@ -354,7 +352,7 @@ void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray trackids.push_back(trid); } EdbTrackP *t_chosen = 0; - float r2, dz, imp, imp2, r2max, dzmax, impmax, impmax2; r2max=dzmax=impmax=impmax2=1e9f; + float r2, dz, r2max, dzmax; r2max=dzmax=1e9f; int founds=0; for(int it=0; itSetFlag(999999); - if( r2 < r2max ) { r2max=r2; dzmax=dz; impmax=imp; impmax2=imp2; t_chosen=t; } + if( r2 < r2max ) { r2max=r2; dzmax=dz; t_chosen=t; } // v_vtx.AddTrack(t); founds++; - outTree->Fill(0, 1, v->ID(), t->ID() ,t->N(), t->Npl(), t->TX(), t->TY(), t->GetSegmentFirst()->Plate(), t->GetSegmentLast()->Plate(), r2, dz, imp, imp2); + outTree->Fill(0, 1, v->ID(), t->ID() ,t->N(), t->Npl(), t->TX(), t->TY(), t->GetSegmentFirst()->Plate(), t->GetSegmentLast()->Plate(), r2, dz); } } if (flag1){ v_vtx.AddTrack(t_chosen); - outTree->Fill(1, founds, v->ID(), t_chosen->ID(), t_chosen->N(), t_chosen->Npl(), t_chosen->TX(), t_chosen->TY(), t_chosen->GetSegmentFirst()->Plate(), t_chosen->GetSegmentLast()->Plate(), r2max, dzmax, impmax, impmax2); + outTree->Fill(1, founds, v->ID(), t_chosen->ID(), t_chosen->N(), t_chosen->Npl(), t_chosen->TX(), t_chosen->TY(), t_chosen->GetSegmentFirst()->Plate(), t_chosen->GetSegmentLast()->Plate(), r2max, dzmax); Log(1,"AddCompatibleTracks","Closest track found at r2=%.4f dz=%.2f\n",r2max,dzmax); v_out2.Add(v); } @@ -382,7 +380,7 @@ void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray } } -bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz, float *imp, float *imp2) +bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz) { EdbSegP ss; EdbSegP *tst = t.GetSegmentFirst(); @@ -393,11 +391,8 @@ bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz, f float dy=ss.Y()-v.VY(); *r2 = Sqrt(dx*dx+dy*dy); *dz = ss.DZ(); - *imp = distance(t, v); - *imp2 = distance(ss, v); float r2max = env.GetValue("emvertex.trfit.r2max" , 5. ); float dzmax = env.GetValue("emvertex.trfit.dzmax" , 4000. ); - float impmax = env.GetValue("emvertex.trfit.impmax" , 5. ); if(*r2 Date: Mon, 24 Mar 2025 19:03:23 +0100 Subject: [PATCH 13/19] Rebase to recent updates --- src/appl/emrec/emvertex.cpp | 1 + 1 file changed, 1 insertion(+) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index 042211d5..c7b82f4e 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -377,6 +377,7 @@ void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray else { v_out.Add(v); } + else { v_out2.Add(v);} } } From ae7d16b97b82e1fd53ea5f60b852b11a1678f2e8 Mon Sep 17 00:00:00 2001 From: dcentanni Date: Mon, 24 Mar 2025 14:48:38 +0100 Subject: [PATCH 14/19] Some fixes --- src/appl/emrec/emvertex.cpp | 1 - 1 file changed, 1 deletion(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index c7b82f4e..042211d5 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -377,7 +377,6 @@ void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray else { v_out.Add(v); } - else { v_out2.Add(v);} } } From c5bd77ac8240398d3c4184058f99716f038ef390 Mon Sep 17 00:00:00 2001 From: dcentanni Date: Fri, 28 Mar 2025 16:07:46 +0100 Subject: [PATCH 15/19] Vertex track association refit for emvertex application --- src/appl/emrec/emvertex.cpp | 75 +++++++++++++++++++++++++++++++++++-- 1 file changed, 71 insertions(+), 4 deletions(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index 042211d5..20da7235 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -20,6 +20,9 @@ EdbID idset; EdbPVRec gAli; EdbScanProc gSproc; EdbVertexRec gEVR; +EdbVertexRec rfEVR; + +bool do_vtxrefit = false; void VertexRec(EdbID id, TEnv &cenv); void ReadVertex(EdbID id,TEnv &env); @@ -28,6 +31,9 @@ void SetTracksErrors(TObjArray &tracks, EdbScanCond &cond, float p, float m); void do_vertex(TEnv &env); void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TObjArray &v_out2, TNtuple* outTree); bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz); +void SplitTrack(EdbTrackP *t, EdbTrackP *&t_in, EdbTrackP *&t_out, Int_t zsplit); +void ExecuteVTA(TEnv &env, EdbVertex *vtx, EdbTrackP *track); +int SetSegmentsP(EdbTrackP t, float p) {for(int i=0; iSetP(p); return t.N();} void Display( const char *dsname, EdbVertexRec *evr, TEnv &env ); //---------------------------------------------------------------------------------------- @@ -167,8 +173,7 @@ int main(int argc, char* argv[]) bool do_set = false; bool do_display = false; bool do_read = false; - bool do_vtxrefit = false; - + for(int i=1; i trackids; EdbVertex *v = v_vtx.GetVertex(iv); + Log(1,"AddCompatibleTracks","Looking for parent tracks of vtx: %i\n",v->ID()); for(int i=0; iN(); i++){ EdbTrackP *t = (EdbTrackP*)v->GetTrack(i); @@ -373,11 +379,26 @@ void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray outTree->Fill(1, founds, v->ID(), t_chosen->ID(), t_chosen->N(), t_chosen->Npl(), t_chosen->TX(), t_chosen->TY(), t_chosen->GetSegmentFirst()->Plate(), t_chosen->GetSegmentLast()->Plate(), r2max, dzmax); Log(1,"AddCompatibleTracks","Closest track found at r2=%.4f dz=%.2f\n",r2max,dzmax); v_out2.Add(v); + if (do_vtxrefit) + { + Log(1,"AddCompatibleTracks","Executing VTA on vertex vID=%d",v->ID()); + rfEVR.SetPVRec(&v_vtx); + ExecuteVTA(env, v, t_chosen); + } } - else { + else + { + + if (do_vtxrefit) {rfEVR.AddVertex(v); cout << "Adding vertex " << v->ID() << " with flag " << v->Flag() << endl;} v_out.Add(v); } } + if (do_vtxrefit) + { + TString name; + gSproc.MakeFileName(name,idset,"vtx.refit.root",false); + EdbDataProc::MakeVertexTree(*(rfEVR.eVTX),name.Data()); + } } bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz) @@ -397,6 +418,52 @@ bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz) return false; } +void SplitTrack(EdbTrackP *t, EdbTrackP *&t_in, EdbTrackP *&t_out, Int_t zsplit) +{ + TEnv trenv("trenv"); + trenv.ReadFile("track.rootrc", kEnvLocal); + float tr_pfit = trenv.GetValue("fedra.track.momentum" , 1000); + float tr_mfit = trenv.GetValue("fedra.track.mass" , 0.14); + EdbSegP *sbest = (EdbSegP *) t->GetSegmentWithClosestZ(zsplit, 3000); + if (!sbest) {Log(1, "EdbTrackP::GetSegmentWithClosestZ", "closest segment not found!"); return;} + Int_t cutplate = -1; + if (sbest->Z() < zsplit) cutplate = sbest->Plate(); + else {cutplate = sbest->Plate() - 1;} + for (Int_t iseg = 0; iseg < t->N(); iseg++) + { + EdbSegP *seg = (EdbSegP *) t->GetSegment(iseg); + if (seg->Plate() <= cutplate) t_in->AddSegment(seg); + else {t_out->AddSegment(seg);} + } + SetSegmentsP(t_in, tr_pfit); + t_in->SetM(tr_mfit); + t_in->SetNpl(); + t_in->SetID(t->ID()); + t_in->FitTrackKFS(); + SetSegmentsP(t_out, tr_pfit); + t_out->SetM(tr_mfit); + t_out->SetNpl(); + t_out->SetID(t->ID()*10); + t_out->FitTrackKFS(); +} +void ExecuteVTA(TEnv &env, EdbVertex *vtx, EdbTrackP *track) +{ + EdbVTA *vta = NULL; + // Make a new EdbVertex object in order to not change the original EdbVertex obj + EdbVertex *vtx_new = new EdbVertex(); + vtx_new->SetV(vtx->V()); + for (int t=0;tN();t++){ vta = rfEVR.AddTrack(*vtx_new, (EdbTrackP*)vtx->GetTrack(t), true); } + EdbTrackP *intrack = new EdbTrackP(); + EdbTrackP *outtrack = new EdbTrackP(); + SplitTrack(track, intrack, outtrack, vtx->VZ()); + vta = rfEVR.AddTrack(*vtx_new, intrack, false); + vta = rfEVR.AddTrack(*vtx_new, outtrack, true); + //vtx_new->EstimateVertexFlag(); + vtx_new->SetID(vtx->ID()); + cout << "Estimated vertex flag for vtx_new ID: " << vtx_new->ID() << " " << vtx_new->Flag() << " " << vtx_new->EstimateVertexFlag() << endl; + rfEVR.AddVertex(vtx_new); +} + //----------------------------------------------------------------------------- void SetTracksErrors(TObjArray &tracks, EdbScanCond &cond, float p, float m) { From a241700853016645ddec824c23b4b762f96c1ca8 Mon Sep 17 00:00:00 2001 From: Fabio Date: Fri, 28 Mar 2025 16:47:16 +0100 Subject: [PATCH 16/19] GetEnv out of loop, debug verbosity --- src/appl/emrec/emvertex.cpp | 46 ++++++++++++++++--------------------- 1 file changed, 20 insertions(+), 26 deletions(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index 20da7235..0f48b687 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -29,10 +29,10 @@ void ReadVertex(EdbID id,TEnv &env); void MakeScanCondBT(EdbScanCond &cond, TEnv &env); void SetTracksErrors(TObjArray &tracks, EdbScanCond &cond, float p, float m); void do_vertex(TEnv &env); -void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TObjArray &v_out2, TNtuple* outTree); -bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz); +void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, float r2max, float dzmax, TObjArray &v_out, TObjArray &v_out2, TNtuple* outTree); +bool IsCompatible(EdbVertex &v, EdbTrackP &t, float r2max, float dzmax, float *r2, float *dz); void SplitTrack(EdbTrackP *t, EdbTrackP *&t_in, EdbTrackP *&t_out, Int_t zsplit); -void ExecuteVTA(TEnv &env, EdbVertex *vtx, EdbTrackP *track); +void ExecuteVTA(EdbVertex *vtx, EdbTrackP *track); int SetSegmentsP(EdbTrackP t, float p) {for(int i=0; iSetP(p); return t.N();} void Display( const char *dsname, EdbVertexRec *evr, TEnv &env ); @@ -242,6 +242,9 @@ void ReadVertex(EdbID id, TEnv &env) gEVR.eQualityMode= env.GetValue("emvertex.vtx.QualityMode" , 0); // (0:=Prob/(sigVX^2+sigVY^2); 1:= inverse average track-vertex distance) TCut cutvtx = env.GetValue("emvertex.vtx.cutvtx" , "(flag==0||flag==3)&&n>4"); + float r2max = env.GetValue("emvertex.trfit.r2max" , 5. ); + float dzmax = env.GetValue("emvertex.trfit.dzmax" , 4000. ); + TObjArray v_out; TObjArray v_out2; @@ -258,7 +261,7 @@ void ReadVertex(EdbID id, TEnv &env) vtr->SetScanCond( new EdbScanCond(gCond) ); gSproc.ReadTracksTree( idset,*vtr, cuttr); TNtuple *outTree = new TNtuple("tracks","Tree of matched tracks","chosen:n:vid:tid:nseg:npl:tx:ty:firstp:lastp:r2:dz"); - AddCompatibleTracks(env, *vtr, gAli , v_out, v_out2, outTree); // assign to the vertices of gAli additional tracks from vtr if any + AddCompatibleTracks(*vtr, gAli, r2max, dzmax, v_out, v_out2, outTree); // assign to the vertices of gAli additional tracks from vtr if any EdbDataProc::MakeVertexTree(v_out,"flag0.vtx.root"); EdbDataProc::MakeVertexTree(v_out2,"flag1.vtx.root"); TFile *outFile = new TFile("found_tracks.root","RECREATE"); @@ -340,7 +343,7 @@ void MakeScanCondBT(EdbScanCond &cond, TEnv &env) cond.SetName("SND_basetrack"); } -void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray &v_out, TObjArray &v_out2, TNtuple* outTree) +void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, float r2max, float dzmax, TObjArray &v_out, TObjArray &v_out2, TNtuple* outTree) { int ntr = v_trk.Ntracks(); int nvtx = v_vtx.Nvtx(); @@ -365,7 +368,7 @@ void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray EdbTrackP *t = v_trk.GetTrack(it); int trid = t->ID(); if (std::find(trackids.begin(), trackids.end(), trid)!=trackids.end()) continue; //Maybe here can be changed to EdbVertex::TrackInVertex(EdbTrackP *t) - if( IsCompatible(env, *v,*t, &r2, &dz) ) { + if( IsCompatible(*v, *t, r2max, dzmax, &r2, &dz) ) { flag1 = true; t->SetFlag(999999); if( r2 < r2max ) { r2max=r2; dzmax=dz; t_chosen=t; } @@ -377,13 +380,13 @@ void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray if (flag1){ v_vtx.AddTrack(t_chosen); outTree->Fill(1, founds, v->ID(), t_chosen->ID(), t_chosen->N(), t_chosen->Npl(), t_chosen->TX(), t_chosen->TY(), t_chosen->GetSegmentFirst()->Plate(), t_chosen->GetSegmentLast()->Plate(), r2max, dzmax); - Log(1,"AddCompatibleTracks","Closest track found at r2=%.4f dz=%.2f\n",r2max,dzmax); + Log(2,"AddCompatibleTracks","Closest track found at r2=%.4f dz=%.2f\n",r2max,dzmax); v_out2.Add(v); if (do_vtxrefit) { - Log(1,"AddCompatibleTracks","Executing VTA on vertex vID=%d",v->ID()); + Log(2,"AddCompatibleTracks","Executing VTA on vertex vID=%d",v->ID()); rfEVR.SetPVRec(&v_vtx); - ExecuteVTA(env, v, t_chosen); + ExecuteVTA(v, t_chosen); } } else @@ -401,7 +404,7 @@ void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, TObjArray } } -bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz) +bool IsCompatible(EdbVertex &v, EdbTrackP &t, float r2max, float dzmax, float *r2, float *dz) { EdbSegP ss; EdbSegP *tst = t.GetSegmentFirst(); @@ -412,9 +415,11 @@ bool IsCompatible(TEnv &env, EdbVertex &v, EdbTrackP &t, float *r2, float *dz) float dy=ss.Y()-v.VY(); *r2 = Sqrt(dx*dx+dy*dy); *dz = ss.DZ(); - float r2max = env.GetValue("emvertex.trfit.r2max" , 5. ); - float dzmax = env.GetValue("emvertex.trfit.dzmax" , 4000. ); - if(*r2SetID(t->ID()*10); t_out->FitTrackKFS(); } -void ExecuteVTA(TEnv &env, EdbVertex *vtx, EdbTrackP *track) +void ExecuteVTA(EdbVertex *vtx, EdbTrackP *track) { EdbVTA *vta = NULL; // Make a new EdbVertex object in order to not change the original EdbVertex obj @@ -481,15 +486,4 @@ void SetTracksErrors(TObjArray &tracks, EdbScanCond &cond, float p, float m) } t->FitTrackKFS(); } -} - -float distance(const EdbTrackP& t, const EdbVertex& v) { - float dx = v.VX() - t.X(); - float dy = v.VY() - t.Y(); - float dz = v.VZ() - t.Z(); - float tx = t.TX(); - float ty = t.TY(); - float nom = std::pow((dx*ty - dy*tx),2) + std::pow((dy - dz*ty),2) + std::pow((dz*tx - dx),2); - float denom = std::pow((tx),2) + std::pow((ty),2) + 1.; -return sqrt(nom/denom); -} +} \ No newline at end of file From 7bebd0934972787d680239ea4bf39b037002bbab Mon Sep 17 00:00:00 2001 From: dcentanni Date: Wed, 2 Apr 2025 14:07:39 +0200 Subject: [PATCH 17/19] Clean up, bug fixes --- src/appl/emrec/emvertex.cpp | 49 +++++++++++++++++++++++-------------- 1 file changed, 31 insertions(+), 18 deletions(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index 0f48b687..d19d1db4 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -23,6 +23,10 @@ EdbVertexRec gEVR; EdbVertexRec rfEVR; bool do_vtxrefit = false; +float tr_pfit = 1000; +float tr_mfit = 0.14; +int last_trkID = -1; + void VertexRec(EdbID id, TEnv &cenv); void ReadVertex(EdbID id,TEnv &env); @@ -244,7 +248,11 @@ void ReadVertex(EdbID id, TEnv &env) float r2max = env.GetValue("emvertex.trfit.r2max" , 5. ); float dzmax = env.GetValue("emvertex.trfit.dzmax" , 4000. ); - + TEnv trenv("trenv"); + trenv.ReadFile("track.rootrc", kEnvLocal); + tr_pfit = trenv.GetValue("fedra.track.momentum" , 1000); + tr_mfit = trenv.GetValue("fedra.track.mass" , 0.14); + TObjArray v_out; TObjArray v_out2; @@ -347,13 +355,20 @@ void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, float r2max, float dz { int ntr = v_trk.Ntracks(); int nvtx = v_vtx.Nvtx(); + rfEVR.SetPVRec(&v_vtx); Log(1,"AddCompatibleTracks", "%d tracks, %d vertex", ntr,nvtx ); + if (do_vtxrefit){ + TFile *trk_f = TFile::Open(Form("b%06d.0.0.0.trk.root", idset.eBrick)); + TTree *trackstree = (TTree*) trk_f->Get("tracks"); + last_trkID = trackstree->GetEntries() - 1; + trk_f->Close(); + } for(int iv=0; ivFlag()==1) {rfEVR.AddVertex(v); continue;} bool flag1 = false; std::vector trackids; - EdbVertex *v = v_vtx.GetVertex(iv); - Log(1,"AddCompatibleTracks","Looking for parent tracks of vtx: %i\n",v->ID()); for(int i=0; iN(); i++){ EdbTrackP *t = (EdbTrackP*)v->GetTrack(i); @@ -361,7 +376,9 @@ void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, float r2max, float dz trackids.push_back(trid); } EdbTrackP *t_chosen = 0; - float r2, dz, r2max, dzmax; r2max=dzmax=1e9f; + float r2, dz; + float r2_ini = r2max; + float dz_ini = dzmax; int founds=0; for(int it=0; itSetFlag(999999); - if( r2 < r2max ) { r2max=r2; dzmax=dz; t_chosen=t; } + if( r2 < r2_ini ) { r2_ini=r2; dz_ini=dz; t_chosen=t; } // v_vtx.AddTrack(t); founds++; outTree->Fill(0, 1, v->ID(), t->ID() ,t->N(), t->Npl(), t->TX(), t->TY(), t->GetSegmentFirst()->Plate(), t->GetSegmentLast()->Plate(), r2, dz); @@ -380,21 +397,20 @@ void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, float r2max, float dz if (flag1){ v_vtx.AddTrack(t_chosen); outTree->Fill(1, founds, v->ID(), t_chosen->ID(), t_chosen->N(), t_chosen->Npl(), t_chosen->TX(), t_chosen->TY(), t_chosen->GetSegmentFirst()->Plate(), t_chosen->GetSegmentLast()->Plate(), r2max, dzmax); - Log(2,"AddCompatibleTracks","Closest track found at r2=%.4f dz=%.2f\n",r2max,dzmax); + Log(2,"AddCompatibleTracks","Closest track found at r2=%.4f dz=%.2f\n",r2_ini,dz_ini); v_out2.Add(v); if (do_vtxrefit) { Log(2,"AddCompatibleTracks","Executing VTA on vertex vID=%d",v->ID()); - rfEVR.SetPVRec(&v_vtx); ExecuteVTA(v, t_chosen); } } else { - if (do_vtxrefit) {rfEVR.AddVertex(v); cout << "Adding vertex " << v->ID() << " with flag " << v->Flag() << endl;} v_out.Add(v); } + //r2max = r2_ini; dzmax = dz_ini; } if (do_vtxrefit) { @@ -417,7 +433,7 @@ bool IsCompatible(EdbVertex &v, EdbTrackP &t, float r2max, float dzmax, float *r *dz = ss.DZ(); if(*r2GetSegmentWithClosestZ(zsplit, 3000); if (!sbest) {Log(1, "EdbTrackP::GetSegmentWithClosestZ", "closest segment not found!"); return;} Int_t cutplate = -1; @@ -443,13 +455,14 @@ void SplitTrack(EdbTrackP *t, EdbTrackP *&t_in, EdbTrackP *&t_out, Int_t zsplit) SetSegmentsP(t_in, tr_pfit); t_in->SetM(tr_mfit); t_in->SetNpl(); - t_in->SetID(t->ID()); + t_in->SetID(last_trkID+1); t_in->FitTrackKFS(); SetSegmentsP(t_out, tr_pfit); t_out->SetM(tr_mfit); t_out->SetNpl(); - t_out->SetID(t->ID()*10); + t_out->SetID(last_trkID+2); t_out->FitTrackKFS(); + last_trkID+=2; } void ExecuteVTA(EdbVertex *vtx, EdbTrackP *track) { @@ -463,10 +476,10 @@ void ExecuteVTA(EdbVertex *vtx, EdbTrackP *track) SplitTrack(track, intrack, outtrack, vtx->VZ()); vta = rfEVR.AddTrack(*vtx_new, intrack, false); vta = rfEVR.AddTrack(*vtx_new, outtrack, true); - //vtx_new->EstimateVertexFlag(); + vtx_new->SetFlag(1); vtx_new->SetID(vtx->ID()); - cout << "Estimated vertex flag for vtx_new ID: " << vtx_new->ID() << " " << vtx_new->Flag() << " " << vtx_new->EstimateVertexFlag() << endl; rfEVR.AddVertex(vtx_new); + Log(2, "SplitTrack", "Track, attached to vertex vID=%d, trid=%d splitted in trid=%d (in) and in trid=%d (out)", vtx->ID(), track->ID(), intrack->ID(), outtrack->ID()); } //----------------------------------------------------------------------------- @@ -486,4 +499,4 @@ void SetTracksErrors(TObjArray &tracks, EdbScanCond &cond, float p, float m) } t->FitTrackKFS(); } -} \ No newline at end of file +} From de82cb312a1b305ca5c9aa7c3b80e3722782e74b Mon Sep 17 00:00:00 2001 From: dcentanni Date: Tue, 6 May 2025 16:04:39 +0200 Subject: [PATCH 18/19] Many bug fixes, working version. Tested on data --- src/appl/emrec/emvertex.cpp | 117 +++++++++++++++++++++++++++--------- 1 file changed, 88 insertions(+), 29 deletions(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index d19d1db4..1933e1eb 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -11,6 +11,8 @@ #include "EdbCombGen.h" #include "EdbVertexComb.h" +#include + using namespace std; using namespace TMath; @@ -24,7 +26,7 @@ EdbVertexRec rfEVR; bool do_vtxrefit = false; float tr_pfit = 1000; -float tr_mfit = 0.14; +float tr_mfit = 0.1390; int last_trkID = -1; @@ -33,7 +35,7 @@ void ReadVertex(EdbID id,TEnv &env); void MakeScanCondBT(EdbScanCond &cond, TEnv &env); void SetTracksErrors(TObjArray &tracks, EdbScanCond &cond, float p, float m); void do_vertex(TEnv &env); -void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, float r2max, float dzmax, TObjArray &v_out, TObjArray &v_out2, TNtuple* outTree); +void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, float r2max, float dzmax, TObjArray &v_out, TObjArray &v_out2, TNtuple* outTree); bool IsCompatible(EdbVertex &v, EdbTrackP &t, float r2max, float dzmax, float *r2, float *dz); void SplitTrack(EdbTrackP *t, EdbTrackP *&t_in, EdbTrackP *&t_out, Int_t zsplit); void ExecuteVTA(EdbVertex *vtx, EdbTrackP *track); @@ -251,7 +253,7 @@ void ReadVertex(EdbID id, TEnv &env) TEnv trenv("trenv"); trenv.ReadFile("track.rootrc", kEnvLocal); tr_pfit = trenv.GetValue("fedra.track.momentum" , 1000); - tr_mfit = trenv.GetValue("fedra.track.mass" , 0.14); + tr_mfit = trenv.GetValue("fedra.track.mass" , 0.1390); TObjArray v_out; TObjArray v_out2; @@ -269,7 +271,7 @@ void ReadVertex(EdbID id, TEnv &env) vtr->SetScanCond( new EdbScanCond(gCond) ); gSproc.ReadTracksTree( idset,*vtr, cuttr); TNtuple *outTree = new TNtuple("tracks","Tree of matched tracks","chosen:n:vid:tid:nseg:npl:tx:ty:firstp:lastp:r2:dz"); - AddCompatibleTracks(*vtr, gAli, r2max, dzmax, v_out, v_out2, outTree); // assign to the vertices of gAli additional tracks from vtr if any + AddCompatibleTracks(env, *vtr, gAli, r2max, dzmax, v_out, v_out2, outTree); // assign to the vertices of gAli additional tracks from vtr if any EdbDataProc::MakeVertexTree(v_out,"flag0.vtx.root"); EdbDataProc::MakeVertexTree(v_out2,"flag1.vtx.root"); TFile *outFile = new TFile("found_tracks.root","RECREATE"); @@ -351,11 +353,18 @@ void MakeScanCondBT(EdbScanCond &cond, TEnv &env) cond.SetName("SND_basetrack"); } -void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, float r2max, float dzmax, TObjArray &v_out, TObjArray &v_out2, TNtuple* outTree) +void AddCompatibleTracks(TEnv &env, EdbPVRec &v_trk, EdbPVRec &v_vtx, float r2max, float dzmax, TObjArray &v_out, TObjArray &v_out2, TNtuple* outTree) { int ntr = v_trk.Ntracks(); int nvtx = v_vtx.Nvtx(); - rfEVR.SetPVRec(&v_vtx); + rfEVR.eEdbTracks = gAli.eTracks; + rfEVR.SetPVRec(&gAli); + rfEVR.eDZmax = env.GetValue("emvertex.vtx.DZmax" , 3000.); + rfEVR.eProbMin = env.GetValue("emvertex.vtx.ProbMinV" , 0.001); + rfEVR.eImpMax = env.GetValue("emvertex.vtx.ImpMax" , 10.); + rfEVR.eUseMom = env.GetValue("emvertex.vtx.UseMom" , false); + rfEVR.eUseSegPar = env.GetValue("emvertex.vtx.UseSegPar" , false); + rfEVR.eQualityMode= env.GetValue("emvertex.vtx.QualityMode" , 0); // (0:=Prob/(sigVX^2+sigVY^2); 1:= inverse average track-vertex distance) Log(1,"AddCompatibleTracks", "%d tracks, %d vertex", ntr,nvtx ); if (do_vtxrefit){ TFile *trk_f = TFile::Open(Form("b%06d.0.0.0.trk.root", idset.eBrick)); @@ -363,13 +372,14 @@ void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, float r2max, float dz last_trkID = trackstree->GetEntries() - 1; trk_f->Close(); } + Log(3, "AddCompatibleTracks", "Last track ID from trk file: %d", last_trkID); for(int iv=0; ivFlag()==1) {rfEVR.AddVertex(v); continue;} bool flag1 = false; std::vector trackids; - Log(1,"AddCompatibleTracks","Looking for parent tracks of vtx: %i\n",v->ID()); + Log(1,"\nAddCompatibleTracks","Looking for parent tracks of vtx: %i",v->ID()); for(int i=0; iN(); i++){ EdbTrackP *t = (EdbTrackP*)v->GetTrack(i); int trid = t->ID(); @@ -383,6 +393,7 @@ void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, float r2max, float dz for(int it=0; itID(); if (std::find(trackids.begin(), trackids.end(), trid)!=trackids.end()) continue; //Maybe here can be changed to EdbVertex::TrackInVertex(EdbTrackP *t) if( IsCompatible(*v, *t, r2max, dzmax, &r2, &dz) ) { @@ -395,9 +406,9 @@ void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, float r2max, float dz } } if (flag1){ - v_vtx.AddTrack(t_chosen); + //v_vtx.AddTrack(t_chosen); outTree->Fill(1, founds, v->ID(), t_chosen->ID(), t_chosen->N(), t_chosen->Npl(), t_chosen->TX(), t_chosen->TY(), t_chosen->GetSegmentFirst()->Plate(), t_chosen->GetSegmentLast()->Plate(), r2max, dzmax); - Log(2,"AddCompatibleTracks","Closest track found at r2=%.4f dz=%.2f\n",r2_ini,dz_ini); + Log(2,"AddCompatibleTracks","Closest track found at r2=%.4f dz=%.2f",r2_ini,dz_ini); v_out2.Add(v); if (do_vtxrefit) { @@ -407,7 +418,7 @@ void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, float r2max, float dz } else { - if (do_vtxrefit) {rfEVR.AddVertex(v); cout << "Adding vertex " << v->ID() << " with flag " << v->Flag() << endl;} + if (do_vtxrefit) {rfEVR.AddVertex(v); Log(2, "AddCompatibleTracks", "Adding vertex %d with flag %d", v->ID(), v->Flag());} v_out.Add(v); } //r2max = r2_ini; dzmax = dz_ini; @@ -417,7 +428,7 @@ void AddCompatibleTracks(EdbPVRec &v_trk, EdbPVRec &v_vtx, float r2max, float dz TString name; gSproc.MakeFileName(name,idset,"vtx.refit.root",false); EdbDataProc::MakeVertexTree(*(rfEVR.eVTX),name.Data()); - } +} } bool IsCompatible(EdbVertex &v, EdbTrackP &t, float r2max, float dzmax, float *r2, float *dz) @@ -430,18 +441,18 @@ bool IsCompatible(EdbVertex &v, EdbTrackP &t, float r2max, float dzmax, float *r float dx=ss.X()-v.VX(); float dy=ss.Y()-v.VY(); *r2 = Sqrt(dx*dx+dy*dy); - *dz = ss.DZ(); + *dz = Abs(ss.DZ()); - if(*r2GetSegmentWithClosestZ(zsplit, 3000); + EdbSegP *sbest = (EdbSegP *) t->GetSegmentWithClosestZ(zsplit, 5000); if (!sbest) {Log(1, "EdbTrackP::GetSegmentWithClosestZ", "closest segment not found!"); return;} Int_t cutplate = -1; if (sbest->Z() < zsplit) cutplate = sbest->Plate(); @@ -452,34 +463,82 @@ void SplitTrack(EdbTrackP *t, EdbTrackP *&t_in, EdbTrackP *&t_out, Int_t zsplit) if (seg->Plate() <= cutplate) t_in->AddSegment(seg); else {t_out->AddSegment(seg);} } - SetSegmentsP(t_in, tr_pfit); + Log(3, "SplitTrack", "Track found is %d, print follows", t->ID()); + if (gEDBDEBUGLEVEL == 3) t->PrintNice(); + //SetSegmentsP(t_in, tr_pfit); + t_in->SetP(tr_pfit); t_in->SetM(tr_mfit); - t_in->SetNpl(); + t_in->SetCounters(); + t_in->SetMC(t->MCEvt(), t->MCTrack()); t_in->SetID(last_trkID+1); - t_in->FitTrackKFS(); - SetSegmentsP(t_out, tr_pfit); + t_in->SetTrack(last_trkID+1); + t_in->SetSegmentsTrack(last_trkID+1); + t_in->SetFlag(999999); + t_in->FitTrackKFS(false, 3504); // using segments with min Z and W radlen + if (t_out->N() != 0) + { + //SetSegmentsP(t_out, tr_pfit); + t_out->SetP(tr_pfit); t_out->SetM(tr_mfit); - t_out->SetNpl(); + t_out->SetCounters(); + t_out->SetMC(t->MCEvt(), t->MCTrack()); t_out->SetID(last_trkID+2); - t_out->FitTrackKFS(); + t_out->SetTrack(last_trkID+2); + t_out->SetSegmentsTrack(last_trkID+2); + t_out->SetFlag(999999); + t_out->FitTrackKFS(false, 3504); last_trkID+=2; + } + else {last_trkID+=1;Log(1, "SplitTrack", "Out-track has 0 segments!");} + Log(3, "SplitTrack", "Track %d is splitted in %d and %d", t->ID(), t_in->ID(), t_out->ID()); + if (gEDBDEBUGLEVEL >= 3) {t_in->PrintNice();t_out->PrintNice();} } void ExecuteVTA(EdbVertex *vtx, EdbTrackP *track) { - EdbVTA *vta = NULL; // Make a new EdbVertex object in order to not change the original EdbVertex obj EdbVertex *vtx_new = new EdbVertex(); - vtx_new->SetV(vtx->V()); - for (int t=0;tN();t++){ vta = rfEVR.AddTrack(*vtx_new, (EdbTrackP*)vtx->GetTrack(t), true); } + Log(2, "ExecuteVTA", "Multiplicity of original vertex is %d", vtx->N()); + //vtx_new->SetV(vtx->V()); // this was causing a bug when accessing vtx->VZ() + for (int t=0;tN();t++) + { + EdbVTA *vta = nullptr; + vta = rfEVR.AddTrack(*vtx_new, (EdbTrackP*)vtx->GetTrack(t), true); + } EdbTrackP *intrack = new EdbTrackP(); EdbTrackP *outtrack = new EdbTrackP(); + if (vtx->N() != vtx_new->N()) + { + Log(2, "ExecuteVTA", "Multiplicities between vtx and vtx_new mismatch before VTA, printing them"); + if( gEDBDEBUGLEVEL >= 2) {vtx->Print();vtx_new->Print();} + } + Log(2, "ExecuteVTA", "Attaching track %d to vertex %d at z=%f", track->ID(), vtx->ID(), vtx->VZ()); SplitTrack(track, intrack, outtrack, vtx->VZ()); - vta = rfEVR.AddTrack(*vtx_new, intrack, false); - vta = rfEVR.AddTrack(*vtx_new, outtrack, true); - vtx_new->SetFlag(1); + EdbVTA *vta_in = new EdbVTA(intrack, vtx_new); + vta_in->SetFlag(2); + vtx_new->AddVTA(vta_in); + vta_in->SetZpos(0); + intrack->AddVTA(vta_in); + EdbVTA *vta_out; + if (outtrack->N()){ + vta_out = new EdbVTA(outtrack, vtx_new); + vta_out->SetFlag(2); + vtx_new->AddVTA(vta_out); + vta_out->SetZpos(1); + outtrack->AddVTA(vta_out); + } + else {Log(1, "ExecuteVTA", "WARNING: Out-track invalid, not attaching it to the vertex");} + vtx_new->EstimateVertexFlag(); + if (vtx_new->Flag() != 1) {Log(2, "ExecuteVTA", "Warning, vtx new %d has flag %d", vtx->ID(), vtx_new->Flag());if (gEDBDEBUGLEVEL == 2) {vtx_new->Print();}} vtx_new->SetID(vtx->ID()); - rfEVR.AddVertex(vtx_new); - Log(2, "SplitTrack", "Track, attached to vertex vID=%d, trid=%d splitted in trid=%d (in) and in trid=%d (out)", vtx->ID(), track->ID(), intrack->ID(), outtrack->ID()); + if (rfEVR.MakeV(*vtx_new) ) {rfEVR.AddVertex(vtx_new);} + else {Log(1, "ExecuteVTA", "VERTEX %d NOT ADDED TO THE VERTEXREC", vtx_new->ID());} + Log(2, "ExecuteVTA", "New vertex created vID=%d, prob is %.3f, the original one was %.3f", vtx_new->ID(), vtx_new->V()->prob(), vtx->V()->prob()); + if (outtrack->N()) {Log(2, "SplitTrack", "Track, attached to vertex vID=%d, trid=%d splitted in trid=%d (in) and in trid=%d (out)", vtx->ID(), track->ID(), intrack->ID(), outtrack->ID());} + else {Log(2, "ExecuteVTA", "Track, attached to vertex vID=%d, trid=%d splitted in trid=%d", vtx->ID(), track->ID(), intrack->ID());} + rfEVR.ePVR->AddTrack(intrack); + if (outtrack->N()) rfEVR.ePVR->AddTrack(outtrack); + SafeDelete(vta_in); + //SafeDelete(vta_out); } //----------------------------------------------------------------------------- From d0a37c8ce26b80431c4a64515b0c5e22c66b81ce Mon Sep 17 00:00:00 2001 From: dcentanni Date: Thu, 8 May 2025 10:47:43 +0200 Subject: [PATCH 19/19] In-track in ExecuteVTA was not saved, fixed --- src/appl/emrec/emvertex.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/appl/emrec/emvertex.cpp b/src/appl/emrec/emvertex.cpp index 1933e1eb..e088a890 100644 --- a/src/appl/emrec/emvertex.cpp +++ b/src/appl/emrec/emvertex.cpp @@ -537,7 +537,7 @@ void ExecuteVTA(EdbVertex *vtx, EdbTrackP *track) else {Log(2, "ExecuteVTA", "Track, attached to vertex vID=%d, trid=%d splitted in trid=%d", vtx->ID(), track->ID(), intrack->ID());} rfEVR.ePVR->AddTrack(intrack); if (outtrack->N()) rfEVR.ePVR->AddTrack(outtrack); - SafeDelete(vta_in); + //SafeDelete(vta_in); //SafeDelete(vta_out); }