Investigating the fit difference between DaVinci versions and reconstruction script versions

[yipengsun]

We've see big discrepancies in various fit variables between the following two combo:

1. DaVinci/v36r1p2 + Phoebe's script
2. DaVinci/v42r8p1 + Yipeng's script

To investigate further, we plan to run:

1. Phoebe's original script
2. Yipeng's script around 2019-09

on DaVinci/v36r1p2 with identical inputs to see the differences.

[yipengsun]

So far, I'm able to:

1. Run Phoebe's script with DV36
2. Run my script with DV36
3. Run my script with DV45

I'm scratching my head trying to make Phoebe's script run for DV45, but I got the following error message, which doesn't make too much sense to me:

SeqYMaker            INFO Member list: LoKi::HDRFilter/TriggeredD0, LoKi::HDRFilter/StrippedBCands, LoKi::VoidFilter/SELECT:/Event/Semileptonic/Phys/b2D0MuXB2DMuNuForTauMuLine, FilterDesktop/smuTISfilterSel, FilterInTrees/chargedK, FilterInTrees/chargedPi, CombineParticles/SelMyD0, LoKi::VoidFilter/SELECT:Phys/StdAllLoosePions, LoKi::VoidFilter/SELECT:Phys/StdNoPIDsUpPions, CombineParticles/SelMyDst, FilterInTrees/mulist, CombineParticles/SelMyBd, FitDecayTrees/YDTFSel
TriggeredD0         ERROR LoKi::HDRFilter:: Inconsistent setting of code&location       StatusCode=FAILURE
TriggeredD0       WARNING LoKi::HDRFilter:: Inconsistent setting of name&location       StatusCode=FAILURE
TriggeredD0       WARNING LoKi::HDRFilter:: Inconsistent setting of name&code&location  StatusCode=FAILURE
SeqYMaker2           INFO Member list: LoKi::HDRFilter/TriggeredD0, LoKi::HDRFilter/StrippedBCands, LoKi::VoidFilter/SELECT:/Event/Semileptonic/Phys/b2D0MuXB2DMuNuForTauMuLine, FilterDesktop/smuTISfilterSel, FilterInTrees/chargedK, FilterInTrees/chargedPi, CombineParticles/SelMyD0, LoKi::VoidFilter/SELECT:Phys/StdAllLoosePions, LoKi::VoidFilter/SELECT:Phys/StdNoPIDsUpPions, CombineParticles/SelMyDst, FilterInTrees/mulist, CombineParticles/SelMyWSBd, FitDecayTrees/YDTFSelWS
SeqYMaker3           INFO Member list: LoKi::HDRFilter/TriggeredD0, LoKi::HDRFilter/StrippedBCands, LoKi::VoidFilter/SELECT:/Event/Semileptonic/Phys/b2D0MuXB2DMuNuForTauMuLine, FilterDesktop/smuTISfilterSel, FilterInTrees/chargedK, FilterInTrees/chargedPi, CombineParticles/SelMyD0, LoKi::VoidFilter/SELECT:Phys/StdAllLoosePions, LoKi::VoidFilter/SELECT:Phys/StdNoPIDsUpPions, CombineParticles/SelMyDst2, FilterInTrees/mulist, CombineParticles/SelMyWS2Bd, FitDecayTrees/YDTFSelWS2
StdNoPIDsVeloPions   INFO Particle/AntiParticle to be created : 'pi+'/'pi-'
YCands               INFO Will look for  ( ( B~0  -> ^( ( D*(2010)+  -> ^( ( D0  -> ^( K- ) ^( pi+ )) ) ^( pi+ )) ) ^( mu- )) || ( B0  -> ^( ( D*(2010)-  -> ^( ( D~0  -> ^( K+ ) ^( pi- )) ) ^( pi- )) ) ^( mu+ )) )
YCands               INFO Explicit sub decays :7
YCands.D0            INFO YCands.D0 : ( ( B0  -> ( D*(2010)-  -> ^( ( D~0  -> K+  pi- ) ) pi- )  mu+ ) || ( B~0  -> ( D*(2010)+  -> ^( ( D0  -> K-  pi+ ) ) pi+ )  mu- ) )
YCands.Dst_2010...   INFO YCands.Dst_2010_minus : ( ( B0  -> ^( ( D*(2010)-  -> ( D~0  -> K+  pi- )  pi- ) ) mu+ ) || ( B~0  -> ^( ( D*(2010)+  -> ( D0  -> K-  pi+ )  pi+ ) ) mu- ) )
YCands.Kplus         INFO YCands.Kplus : ( ( B0  -> ( D*(2010)-  -> ( D~0  -> ^( K+ ) pi- )  pi- )  mu+ ) || ( B~0  -> ( D*(2010)+  -> ( D0  -> ^( K- ) pi+ )  pi+ )  mu- ) )
YCands.Y             INFO YCands.Y : ^( ( ( B0  -> ( D*(2010)-  -> ( D~0  -> K+  pi- )  pi- )  mu+ ) || ( B~0  -> ( D*(2010)+  -> ( D0  -> K-  pi+ )  pi+ )  mu- ) ) )
YCands.muplus        INFO YCands.muplus : ( ( B0  -> ( D*(2010)-  -> ( D~0  -> K+  pi- )  pi- )  ^( mu+ )) || ( B~0  -> ( D*(2010)+  -> ( D0  -> K-  pi+ )  pi+ )  ^( mu- )) )
YCands.piminus       INFO YCands.piminus : ( ( B0  -> ( D*(2010)-  -> ( D~0  -> K+  pi- )  ^( pi- ))  mu+ ) || ( B~0  -> ( D*(2010)+  -> ( D0  -> K-  pi+ )  ^( pi+ ))  mu- ) )
YCands.piminus0      INFO YCands.piminus0 : ( ( B0  -> ( D*(2010)-  -> ( D~0  -> K+  ^( pi- ))  pi- )  mu+ ) || ( B~0  -> ( D*(2010)+  -> ( D0  -> K-  ^( pi+ ))  pi+ )  mu- ) )
YCandsWS             INFO Will look for  ( ( B~0  -> ^( ( D*(2010)+  -> ^( ( D0  -> ^( K- ) ^( pi+ )) ) ^( pi+ )) ) ^( mu+ )) || ( B0  -> ^( ( D*(2010)-  -> ^( ( D~0  -> ^( K+ ) ^( pi- )) ) ^( pi- )) ) ^( mu- )) )
YCandsWS             INFO Explicit sub decays :7
YCandsWS.D0          INFO YCandsWS.D0 : ( ( B0  -> ( D*(2010)-  -> ^( ( D~0  -> K+  pi- ) ) pi- )  mu- ) || ( B~0  -> ( D*(2010)+  -> ^( ( D0  -> K-  pi+ ) ) pi+ )  mu+ ) )
YCandsWS.Dst_20...   INFO YCandsWS.Dst_2010_minus : ( ( B0  -> ^( ( D*(2010)-  -> ( D~0  -> K+  pi- )  pi- ) ) mu- ) || ( B~0  -> ^( ( D*(2010)+  -> ( D0  -> K-  pi+ )  pi+ ) ) mu+ ) )
YCandsWS.Kplus       INFO YCandsWS.Kplus : ( ( B0  -> ( D*(2010)-  -> ( D~0  -> ^( K+ ) pi- )  pi- )  mu- ) || ( B~0  -> ( D*(2010)+  -> ( D0  -> ^( K- ) pi+ )  pi+ )  mu+ ) )
YCandsWS.Y           INFO YCandsWS.Y : ^( ( ( B0  -> ( D*(2010)-  -> ( D~0  -> K+  pi- )  pi- )  mu- ) || ( B~0  -> ( D*(2010)+  -> ( D0  -> K-  pi+ )  pi+ )  mu+ ) ) )
YCandsWS.muplus      INFO YCandsWS.muplus : ( ( B0  -> ( D*(2010)-  -> ( D~0  -> K+  pi- )  pi- )  ^( mu- )) || ( B~0  -> ( D*(2010)+  -> ( D0  -> K-  pi+ )  pi+ )  ^( mu+ )) )
YCandsWS.piminus     INFO YCandsWS.piminus : ( ( B0  -> ( D*(2010)-  -> ( D~0  -> K+  pi- )  ^( pi- ))  mu- ) || ( B~0  -> ( D*(2010)+  -> ( D0  -> K-  pi+ )  ^( pi+ ))  mu+ ) )
YCandsWS.piminus0    INFO YCandsWS.piminus0 : ( ( B0  -> ( D*(2010)-  -> ( D~0  -> K+  ^( pi- ))  pi- )  mu- ) || ( B~0  -> ( D*(2010)+  -> ( D0  -> K-  ^( pi+ ))  pi+ )  mu+ ) )
YCandsWS2            INFO Will look for  ( ( B~0  -> ^( ( D*(2010)-  -> ^( ( D0  -> ^( K- ) ^( pi+ )) ) ^( pi- )) ) ^( mu- )) || ( B0  -> ^( ( D*(2010)+  -> ^( ( D~0  -> ^( K+ ) ^( pi- )) ) ^( pi+ )) ) ^( mu+ )) )
YCandsWS2            INFO Explicit sub decays :7
YCandsWS2.D0         INFO YCandsWS2.D0 : ( ( B0  -> ( D*(2010)+  -> ^( ( D~0  -> K+  pi- ) ) pi+ )  mu+ ) || ( B~0  -> ( D*(2010)-  -> ^( ( D0  -> K-  pi+ ) ) pi- )  mu- ) )
YCandsWS2.Dst_2...   INFO YCandsWS2.Dst_2010_minus : ( ( B0  -> ^( ( D*(2010)+  -> ( D~0  -> K+  pi- )  pi+ ) ) mu+ ) || ( B~0  -> ^( ( D*(2010)-  -> ( D0  -> K-  pi+ )  pi- ) ) mu- ) )
YCandsWS2.Kplus      INFO YCandsWS2.Kplus : ( ( B0  -> ( D*(2010)+  -> ( D~0  -> ^( K+ ) pi- )  pi+ )  mu+ ) || ( B~0  -> ( D*(2010)-  -> ( D0  -> ^( K- ) pi+ )  pi- )  mu- ) )
YCandsWS2.Y          INFO YCandsWS2.Y : ^( ( ( B0  -> ( D*(2010)+  -> ( D~0  -> K+  pi- )  pi+ )  mu+ ) || ( B~0  -> ( D*(2010)-  -> ( D0  -> K-  pi+ )  pi- )  mu- ) ) )
YCandsWS2.muplus     INFO YCandsWS2.muplus : ( ( B0  -> ( D*(2010)+  -> ( D~0  -> K+  pi- )  pi+ )  ^( mu+ )) || ( B~0  -> ( D*(2010)-  -> ( D0  -> K-  pi+ )  pi- )  ^( mu- )) )
YCandsWS2.piminus    INFO YCandsWS2.piminus : ( ( B0  -> ( D*(2010)+  -> ( D~0  -> K+  pi- )  ^( pi+ ))  mu+ ) || ( B~0  -> ( D*(2010)-  -> ( D0  -> K-  pi+ )  ^( pi- ))  mu- ) )
YCandsWS2.piminus0   INFO YCandsWS2.piminus0 : ( ( B0  -> ( D*(2010)+  -> ( D~0  -> K+  ^( pi- ))  pi+ )  mu+ ) || ( B~0  -> ( D*(2010)-  -> ( D0  -> K-  ^( pi+ ))  pi- )  mu- ) )
MonitoringSequence   INFO Member list:
LumiSeq              INFO OR Member list: EventAccounting/EventAccount
EventPersistenc...   INFO Added successfully Conversion service:RootCnvSvc
EventSelector        INFO Stream:EventSelector.DataStreamTool_1 Def:DATAFILE='./data/data-2012-md/00041836_00006100_1.semileptonic.dst' SVC='Gaudi::RootEvtSelector' OPT='READ'
ApplicationMgr       INFO Application Manager Initialized successfully
ApplicationMgr       INFO Application Manager Started successfully
FileRecordPersi...   INFO Added successfully Conversion service:FileRecordCnvSvc
EventSelector     SUCCESS Reading Event record 1. Record number within stream 1: 1
StateScale           INFO  Condition: /MomentumScale     Valid;   Validity: Sat Dec 31 23:00:00 2011 -> Sat Dec  1 00:00:00 2012
MagneticFieldSvc     INFO Map scaled by factor 0.999998 with polarity internally used: -1 signed relative current: -0.999998
TriggeredD0         FATAL LoKi::HDRFilter:: Exception throw: get():: No valid data at 'Hlt/DecReports' StatusCode=FAILURE
TriggeredD0         FATAL  Exception with tag= is caught
TriggeredD0         ERROR    TriggeredD0:: get():: No valid data at 'Hlt/DecReports'     StatusCode=FAILURE
SeqYMaker           FATAL  Exception with tag= is caught
SeqYMaker           ERROR    TriggeredD0:: get():: No valid data at 'Hlt/DecReports'     StatusCode=FAILURE
DaVinciUserSequ...  FATAL  Exception with tag= is caught
DaVinciUserSequ...  ERROR    TriggeredD0:: get():: No valid data at 'Hlt/DecReports'     StatusCode=FAILURE
DaVinciAnalysisSeq  FATAL  Exception with tag= is caught
DaVinciAnalysisSeq  ERROR    TriggeredD0:: get():: No valid data at 'Hlt/DecReports'     StatusCode=FAILURE
FilteredEventSeq    FATAL  Exception with tag= is caught
FilteredEventSeq    ERROR    TriggeredD0:: get():: No valid data at 'Hlt/DecReports'     StatusCode=FAILURE
DaVinciEventSeq     FATAL  Exception with tag= is caught
DaVinciEventSeq     ERROR    TriggeredD0:: get():: No valid data at 'Hlt/DecReports'     StatusCode=FAILURE
EventLoopMgr        FATAL .executeEvent(): Exception with tag= thrown by DaVinciEventSeq
EventLoopMgr        ERROR    TriggeredD0:: get():: No valid data at 'Hlt/DecReports'     StatusCode=FAILURE
EventLoopMgr      WARNING Execution of algorithm DaVinciEventSeq failed
EventLoopMgr        ERROR Error processing event loop.
EventLoopMgr        ERROR Terminating event processing loop due to errors
EventLoopMgr        ERROR Terminating event processing loop due to errors
ApplicationMgr       INFO Application Manager Stopped successfully

As far as I can tell, the definition of TriggeredD0 is the same in both scripts.

[yipengsun]

I can confirm that if I remove TriggeredD0 (named as trigfltr in Phoebe's script), then Phoebe's script can run with DV45. I just can't see how the trigfltr is defined incorrectly!

[yipengsun]

I've finished generating all required ntuples, and pretty much give up on making Phoebe's old script work with DV45 (see 2 posts above for more details).

The info on DaVinci scripts used for this test, and how to run them on lxplus, are documented in README of test-comp-davinci branch.

The generated ntuples can be found in the master branch, inside run1-b2D0MuXB2DMuNuForTauMuLine/samples folder:

Dst--20_08_07--std--data--2012--md--dv36-subset-phoebe.root
Dst--20_08_07--std--data--2012--md--dv36-subset-yipeng.root
Dst--20_08_10--std--data--2012--md--dv45-subset-yipeng.root

Use git annex to get them.

The corresponding DaVinci logs are:

• Dst-20_08_07-std-dv36_phoebe.log
• Dst-20_08_07-std-dv36_yipeng.log
• Dst-20_08_10-std-dv45_yipeng.log

FYI @manuelfs

[manuelfs]

Nice job! With some fortune, by tomorrow you will have found the origin of the difference between our ntuple and Phoebe's.

By the way, why would we want to run Phoebe's in DV45?

[yipengsun]

Just so we can have more comparison, if needed

[yipengsun]

Here's some comparison plots with latest ntuples:

DV36 Phoebe vs. DV36 Yipeng

DV36 Phoebe vs. DV45 Yipeng

In both cases, we cannot see the discrepancies observed in the large ntuple comparison.

[yipengsun]

My plan to move forward:

1. Still working on making the reco_Dst-debug-yipeng.py script produces identical result with reco_Dst-debug-phoebe.py. By doing this, I'll be able to find the subtle differences and have a better understanding of DaVinci

2. Phoebe's large ntuple is still not reproducible.

3. I think we had 2 important questions:

   1. Why do we seem to have one extra constraint when doing the fit?
   2. Why we cannot reconstruct a D0 Mu combo for every event that passes stripping

I don't think we have answer to either. Given that we still cannot see the first problems with our small ntuples, I don't know what to do for the first. For the second, we should keep asking.

[yipengsun]

After removing the upstream pions from Phoebe's, the results are much closer

DV36 Phoebe vs. DV36 Yipeng

DV36 Phoebe vs. DV45 Yipeng

[yipengsun]

I think we FINALLY FOUND THE PROBLEM! Now I only refit D, instead of D Mu combo, in Phoebe's script, and I did another comparison.

DV36 Phoebe, refit D* only vs. DV36 Yipeng

[yipengsun]

I generated component-wise cutflows with latest cocktail ntuples: https://umd-lhcb.github.io/lhcb-ntuples-gen/data/cutflows/Dst_refit_Dst_only_vs_full_refit/

We lose events in signal modes when doing a FULL refit, comparing to refit D* only, but we have more events in normalization and D** modes.

Also, the full refit vs refit D* only comparison plots are consistent with what we saw in data ntuples: