Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
D0TruthSeparator.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file D0TruthSeparator.C
1 #include <cstdlib>
2 #include <iostream>
3 #include <map>
4 #include <string>
5 #include <vector>
6 
7 #include "TChain.h"
8 #include "TFile.h"
9 #include "TTree.h"
10 #include "TString.h"
11 #include "TObjString.h"
12 #include "TSystem.h"
13 #include "TROOT.h"
14 
15 void D0TruthSeparator()
16 {
17  string inDir = "/sphenix/user/rosstom/analysis/HF-Particle/KFParticle_sPHENIX/Run40Acceptance082922/";
18  TString truthFilePath;
19  TString KFPFilePath;
20 
21  Float_t K_mass = .493677;
22  Float_t Pi_mass = .139570;
23 
24  // Initialize new file for saving the separated D0/D0bar information
25  TFile *newfile = new TFile("Run40_D0_Separated_091922.root","recreate");
26  TTree *D0_tree = new TTree("D0_tree","D0_tree");
27  TTree *D0bar_tree = new TTree("D0bar_tree","D0bar_tree");
28  TTree *Background_tree = new TTree("Background_tree","Background_tree");
29 
30  Float_t outTrue_mother_px, outTrue_mother_py, outTrue_mother_pz, outTrue_mother_pE;
31  Float_t outTrue_mother_pT, outTrue_mother_eta;
32  Int_t outTrue_mother_barcode, outTrue_mother_PDG_ID;
33  Float_t outTrue_positive_px, outTrue_positive_py, outTrue_positive_pz, outTrue_positive_pE;
34  Float_t outTrue_negative_px, outTrue_negative_py, outTrue_negative_pz, outTrue_negative_pE;
35  Float_t outTrue_positive_eta, outTrue_negative_eta;
36  Int_t outTrue_positive_PDG_ID, outTrue_negative_PDG_ID;
37  Float_t outKFP_D0_px, outKFP_D0_py, outKFP_D0_pz, outKFP_D0_pE;
38  Float_t outKFP_D0_mass;
39  Float_t outKFP_D0_decayTime, outKFP_D0_decayLength, outKFP_D0_DIRA, outKFP_D0_IP, outKFP_D0_IPchi2;
40  Float_t outKFP_D0_pseudorapidity, outKFP_D0_rapidity;
41  Float_t outKFP_positive_px, outKFP_positive_py, outKFP_positive_pz, outKFP_positive_pE;
42  Float_t outKFP_negative_px, outKFP_negative_py, outKFP_negative_pz, outKFP_negative_pE;
43  Float_t outKFP_positive_p, outKFP_negative_p;
44  Float_t outKFP_positive_pseudorapidity, outKFP_negative_pseudorapidity;
45  Int_t outKFP_positive_charge, outKFP_negative_charge;
46  Float_t outKFP_positive_IP, outKFP_negative_IP;
47  Float_t outKFP_positive_IPchi2, outKFP_negative_IPchi2;
48  Float_t outKFP_track_1_track_2_DCA;
49  Float_t outKFP_KpPm_invm, outKFP_KmPp_invm;
50 
51  D0_tree->Branch("outTrue_mother_px", &outTrue_mother_px,"outTrue_mother_px/F");
52  D0_tree->Branch("outTrue_mother_py", &outTrue_mother_py,"outTrue_mother_py/F");
53  D0_tree->Branch("outTrue_mother_pz", &outTrue_mother_pz,"outTrue_mother_pz/F");
54  D0_tree->Branch("outTrue_mother_pE", &outTrue_mother_pE,"outTrue_mother_pE/F");
55  D0_tree->Branch("outTrue_mother_pT", &outTrue_mother_pT,"outTrue_mother_pT/F");
56  D0_tree->Branch("outTrue_mother_eta", &outTrue_mother_eta,"outTrue_mother_eta/F");
57  D0_tree->Branch("outTrue_mother_barcode", &outTrue_mother_barcode,"outTrue_mother_barcode/I");
58  D0_tree->Branch("outTrue_mother_PDG_ID", &outTrue_mother_PDG_ID,"outTrue_mother_PDG_ID/I");
59  D0_tree->Branch("outTrue_positive_px", &outTrue_positive_px,"outTrue_positive_px/F");
60  D0_tree->Branch("outTrue_positive_py", &outTrue_positive_py,"outTrue_positive_py/F");
61  D0_tree->Branch("outTrue_positive_pz", &outTrue_positive_pz,"outTrue_positive_pz/F");
62  D0_tree->Branch("outTrue_positive_pE", &outTrue_positive_pE,"outTrue_positive_pE/F");
63  D0_tree->Branch("outTrue_negative_px", &outTrue_negative_px,"outTrue_negative_px/F");
64  D0_tree->Branch("outTrue_negative_py", &outTrue_negative_py,"outTrue_negative_py/F");
65  D0_tree->Branch("outTrue_negative_pz", &outTrue_negative_pz,"outTrue_negative_pz/F");
66  D0_tree->Branch("outTrue_negative_pE", &outTrue_negative_pE,"outTrue_negative_pE/F");
67  D0_tree->Branch("outTrue_positive_eta", &outTrue_positive_eta,"outTrue_positive_eta/F");
68  D0_tree->Branch("outTrue_negative_eta", &outTrue_negative_eta,"outTrue_negative_eta/F");
69  D0_tree->Branch("outTrue_positive_PDG_ID", &outTrue_positive_PDG_ID,"outTrue_positive_PDG_ID/I");
70  D0_tree->Branch("outTrue_negative_PDG_ID", &outTrue_negative_PDG_ID,"outTrue_negative_PDG_ID/I");
71  D0_tree->Branch("outKFP_D0_px", &outKFP_D0_px,"outKFP_D0_px/F");
72  D0_tree->Branch("outKFP_D0_py", &outKFP_D0_py,"outKFP_D0_py/F");
73  D0_tree->Branch("outKFP_D0_pz", &outKFP_D0_pz,"outKFP_D0_pz/F");
74  D0_tree->Branch("outKFP_D0_pE", &outKFP_D0_pE,"outKFP_D0_pE/F");
75  D0_tree->Branch("outKFP_D0_mass", &outKFP_D0_mass,"outKFP_D0_mass/F");
76  D0_tree->Branch("outKFP_D0_decayTime", &outKFP_D0_decayTime,"outKFP_D0_decayTime/F");
77  D0_tree->Branch("outKFP_D0_decayLength", &outKFP_D0_decayLength,"outKFP_D0_decayLength/F");
78  D0_tree->Branch("outKFP_D0_DIRA", &outKFP_D0_DIRA,"outKFP_D0_DIRA/F");
79  D0_tree->Branch("outKFP_D0_IP", &outKFP_D0_IP,"outKFP_D0_IP/F");
80  D0_tree->Branch("outKFP_D0_IPchi2", &outKFP_D0_IPchi2,"outKFP_D0_IPchi2/F");
81  D0_tree->Branch("outKFP_D0_pseudorapidity", &outKFP_D0_pseudorapidity,"outKFP_D0_pseudorapidity/F");
82  D0_tree->Branch("outKFP_D0_rapidity", &outKFP_D0_rapidity,"outKFP_D0_rapidity/F");
83  D0_tree->Branch("outKFP_positive_px", &outKFP_positive_px,"outKFP_positive_px/F");
84  D0_tree->Branch("outKFP_positive_py", &outKFP_positive_py,"outKFP_positive_py/F");
85  D0_tree->Branch("outKFP_positive_pz", &outKFP_positive_pz,"outKFP_positive_pz/F");
86  D0_tree->Branch("outKFP_positive_pE", &outKFP_positive_pE,"outKFP_positive_pE/F");
87  D0_tree->Branch("outKFP_negative_px", &outKFP_negative_px,"outKFP_negative_px/F");
88  D0_tree->Branch("outKFP_negative_py", &outKFP_negative_py,"outKFP_negative_py/F");
89  D0_tree->Branch("outKFP_negative_pz", &outKFP_negative_pz,"outKFP_negative_pz/F");
90  D0_tree->Branch("outKFP_negative_pE", &outKFP_negative_pE,"outKFP_negative_pE/F");
91  D0_tree->Branch("outKFP_positive_p", &outKFP_positive_p,"outKFP_positive_p/F");
92  D0_tree->Branch("outKFP_negative_p", &outKFP_negative_p,"outKFP_negative_p/F");
93  D0_tree->Branch("outKFP_positive_pseudorapidity", &outKFP_positive_pseudorapidity,"outKFP_positive_pseudorapidity/F");
94  D0_tree->Branch("outKFP_negative_pseudorapidity", &outKFP_negative_pseudorapidity,"outKFP_negative_pseudorapidity/F");
95  D0_tree->Branch("outKFP_positive_charge", &outKFP_positive_charge,"outKFP_positive_charge/I");
96  D0_tree->Branch("outKFP_negative_charge", &outKFP_negative_charge,"outKFP_negative_charge/I");
97  D0_tree->Branch("outKFP_positive_IP", &outKFP_positive_IP,"outKFP_positive_IP/F");
98  D0_tree->Branch("outKFP_negative_IP", &outKFP_negative_IP,"outKFP_negative_IP/F");
99  D0_tree->Branch("outKFP_positive_IPchi2", &outKFP_positive_IPchi2,"outKFP_positive_IPchi2/F");
100  D0_tree->Branch("outKFP_negative_IPchi2", &outKFP_negative_IPchi2,"outKFP_negative_IPchi2/F");
101  D0_tree->Branch("outKFP_track_1_track_2_DCA", &outKFP_track_1_track_2_DCA,"outKFP_track_1_track_2_DCA/F");
102  D0_tree->Branch("outKFP_KpPm_invm", &outKFP_KpPm_invm,"outKFP_KpPm_invm/F");
103  D0_tree->Branch("outKFP_KmPp_invm", &outKFP_KmPp_invm,"outKFP_KmPp_invm/F");
104 
105  D0bar_tree->Branch("outTrue_mother_px", &outTrue_mother_px,"outTrue_mother_px/F");
106  D0bar_tree->Branch("outTrue_mother_py", &outTrue_mother_py,"outTrue_mother_py/F");
107  D0bar_tree->Branch("outTrue_mother_pz", &outTrue_mother_pz,"outTrue_mother_pz/F");
108  D0bar_tree->Branch("outTrue_mother_pE", &outTrue_mother_pE,"outTrue_mother_pE/F");
109  D0bar_tree->Branch("outTrue_mother_pT", &outTrue_mother_pT,"outTrue_mother_pT/F");
110  D0bar_tree->Branch("outTrue_mother_eta", &outTrue_mother_eta,"outTrue_mother_eta/F");
111  D0bar_tree->Branch("outTrue_mother_barcode", &outTrue_mother_barcode,"outTrue_mother_barcode/I");
112  D0bar_tree->Branch("outTrue_mother_PDG_ID", &outTrue_mother_PDG_ID,"outTrue_mother_PDG_ID/I");
113  D0bar_tree->Branch("outTrue_positive_px", &outTrue_positive_px,"outTrue_positive_px/F");
114  D0bar_tree->Branch("outTrue_positive_py", &outTrue_positive_py,"outTrue_positive_py/F");
115  D0bar_tree->Branch("outTrue_positive_pz", &outTrue_positive_pz,"outTrue_positive_pz/F");
116  D0bar_tree->Branch("outTrue_positive_pE", &outTrue_positive_pE,"outTrue_positive_pE/F");
117  D0bar_tree->Branch("outTrue_negative_px", &outTrue_negative_px,"outTrue_negative_px/F");
118  D0bar_tree->Branch("outTrue_negative_py", &outTrue_negative_py,"outTrue_negative_py/F");
119  D0bar_tree->Branch("outTrue_negative_pz", &outTrue_negative_pz,"outTrue_negative_pz/F");
120  D0bar_tree->Branch("outTrue_negative_pE", &outTrue_negative_pE,"outTrue_negative_pE/F");
121  D0bar_tree->Branch("outTrue_positive_eta", &outTrue_positive_eta,"outTrue_positive_eta/F");
122  D0bar_tree->Branch("outTrue_negative_eta", &outTrue_negative_eta,"outTrue_negative_eta/F");
123  D0bar_tree->Branch("outTrue_positive_PDG_ID", &outTrue_positive_PDG_ID,"outTrue_positive_PDG_ID/I");
124  D0bar_tree->Branch("outTrue_negative_PDG_ID", &outTrue_negative_PDG_ID,"outTrue_negative_PDG_ID/I");
125  D0bar_tree->Branch("outKFP_D0_px", &outKFP_D0_px,"outKFP_D0_px/F");
126  D0bar_tree->Branch("outKFP_D0_py", &outKFP_D0_py,"outKFP_D0_py/F");
127  D0bar_tree->Branch("outKFP_D0_pz", &outKFP_D0_pz,"outKFP_D0_pz/F");
128  D0bar_tree->Branch("outKFP_D0_pE", &outKFP_D0_pE,"outKFP_D0_pE/F");
129  D0bar_tree->Branch("outKFP_D0_mass", &outKFP_D0_mass,"outKFP_D0_mass/F");
130  D0bar_tree->Branch("outKFP_D0_decayTime", &outKFP_D0_decayTime,"outKFP_D0_decayTime/F");
131  D0bar_tree->Branch("outKFP_D0_decayLength", &outKFP_D0_decayLength,"outKFP_D0_decayLength/F");
132  D0bar_tree->Branch("outKFP_D0_DIRA", &outKFP_D0_DIRA,"outKFP_D0_DIRA/F");
133  D0bar_tree->Branch("outKFP_D0_IP", &outKFP_D0_IP,"outKFP_D0_IP/F");
134  D0bar_tree->Branch("outKFP_D0_IPchi2", &outKFP_D0_IPchi2,"outKFP_D0_IPchi2/F");
135  D0bar_tree->Branch("outKFP_D0_pseudorapidity", &outKFP_D0_pseudorapidity,"outKFP_D0_pseudorapidity/F");
136  D0bar_tree->Branch("outKFP_D0_rapidity", &outKFP_D0_rapidity,"outKFP_D0_rapidity/F");
137  D0bar_tree->Branch("outKFP_positive_px", &outKFP_positive_px,"outKFP_positive_px/F");
138  D0bar_tree->Branch("outKFP_positive_py", &outKFP_positive_py,"outKFP_positive_py/F");
139  D0bar_tree->Branch("outKFP_positive_pz", &outKFP_positive_pz,"outKFP_positive_pz/F");
140  D0bar_tree->Branch("outKFP_positive_pE", &outKFP_positive_pE,"outKFP_positive_pE/F");
141  D0bar_tree->Branch("outKFP_negative_px", &outKFP_negative_px,"outKFP_negative_px/F");
142  D0bar_tree->Branch("outKFP_negative_py", &outKFP_negative_py,"outKFP_negative_py/F");
143  D0bar_tree->Branch("outKFP_negative_pz", &outKFP_negative_pz,"outKFP_negative_pz/F");
144  D0bar_tree->Branch("outKFP_negative_pE", &outKFP_negative_pE,"outKFP_negative_pE/F");
145  D0bar_tree->Branch("outKFP_positive_p", &outKFP_positive_p,"outKFP_positive_p/F");
146  D0bar_tree->Branch("outKFP_negative_p", &outKFP_negative_p,"outKFP_negative_p/F");
147  D0bar_tree->Branch("outKFP_positive_pseudorapidity", &outKFP_positive_pseudorapidity,"outKFP_positive_pseudorapidity/F");
148  D0bar_tree->Branch("outKFP_negative_pseudorapidity", &outKFP_negative_pseudorapidity,"outKFP_negative_pseudorapidity/F");
149  D0bar_tree->Branch("outKFP_positive_charge", &outKFP_positive_charge,"outKFP_positive_charge/I");
150  D0bar_tree->Branch("outKFP_negative_charge", &outKFP_negative_charge,"outKFP_negative_charge/I");
151  D0bar_tree->Branch("outKFP_positive_IP", &outKFP_positive_IP,"outKFP_positive_IP/F");
152  D0bar_tree->Branch("outKFP_negative_IP", &outKFP_negative_IP,"outKFP_negative_IP/F");
153  D0bar_tree->Branch("outKFP_positive_IPchi2", &outKFP_positive_IPchi2,"outKFP_positive_IPchi2/F");
154  D0bar_tree->Branch("outKFP_negative_IPchi2", &outKFP_negative_IPchi2,"outKFP_negative_IPchi2/F");
155  D0bar_tree->Branch("outKFP_track_1_track_2_DCA", &outKFP_track_1_track_2_DCA,"outKFP_track_1_track_2_DCA/F");
156  D0bar_tree->Branch("outKFP_KpPm_invm", &outKFP_KpPm_invm,"outKFP_KpPm_invm/F");
157  D0bar_tree->Branch("outKFP_KmPp_invm", &outKFP_KmPp_invm,"outKFP_KmPp_invm/F");
158 
159  //Background_tree->Branch("outTrue_mother_px", &outTrue_mother_px,"outTrue_mother_px/F");
160  //Background_tree->Branch("outTrue_mother_py", &outTrue_mother_py,"outTrue_mother_py/F");
161  //Background_tree->Branch("outTrue_mother_pz", &outTrue_mother_pz,"outTrue_mother_pz/F");
162  //Background_tree->Branch("outTrue_mother_pE", &outTrue_mother_pE,"outTrue_mother_pE/F");
163  //Background_tree->Branch("outTrue_mother_pT", &outTrue_mother_pT,"outTrue_mother_pT/F");
164  //Background_tree->Branch("outTrue_mother_eta", &outTrue_mother_eta,"outTrue_mother_eta/F");
165  //Background_tree->Branch("outTrue_mother_barcode", &outTrue_mother_barcode,"outTrue_mother_barcode/I");
166  //Background_tree->Branch("outTrue_track_1_px", &outTrue_track_1_px,"outTrue_track_1_px/F");
167  //Background_tree->Branch("outTrue_track_1_py", &outTrue_track_1_py,"outTrue_track_1_py/F");
168  //Background_tree->Branch("outTrue_track_1_pz", &outTrue_track_1_pz,"outTrue_track_1_pz/F");
169  //Background_tree->Branch("outTrue_track_1_pE", &outTrue_track_1_pE,"outTrue_track_1_pE/F");
170  //Background_tree->Branch("outTrue_track_2_px", &outTrue_track_2_px,"outTrue_track_2_px/F");
171  //Background_tree->Branch("outTrue_track_2_py", &outTrue_track_2_py,"outTrue_track_2_py/F");
172  //Background_tree->Branch("outTrue_track_2_pz", &outTrue_track_2_pz,"outTrue_track_2_pz/F");
173  //Background_tree->Branch("outTrue_track_2_pE", &outTrue_track_2_pE,"outTrue_track_2_pE/F");
174  //Background_tree->Branch("outTrue_track_1_eta", &outTrue_track_1_eta,"outTrue_track_1_eta/F");
175  //Background_tree->Branch("outTrue_track_2_eta", &outTrue_track_2_eta,"outTrue_track_2_eta/F");
176  //Background_tree->Branch("outTrue_track_1_PDG_ID", &outTrue_track_1_PDG_ID,"outTrue_track_1_PDG_ID/I");
177  //Background_tree->Branch("outTrue_track_2_PDG_ID", &outTrue_track_2_PDG_ID,"outTrue_track_2_PDG_ID/I");
178  Background_tree->Branch("outKFP_D0_px", &outKFP_D0_px,"outKFP_D0_px/F");
179  Background_tree->Branch("outKFP_D0_py", &outKFP_D0_py,"outKFP_D0_py/F");
180  Background_tree->Branch("outKFP_D0_pz", &outKFP_D0_pz,"outKFP_D0_pz/F");
181  Background_tree->Branch("outKFP_D0_pE", &outKFP_D0_pE,"outKFP_D0_pE/F");
182  Background_tree->Branch("outKFP_D0_mass", &outKFP_D0_mass,"outKFP_D0_mass/F");
183  Background_tree->Branch("outKFP_D0_decayTime", &outKFP_D0_decayTime,"outKFP_D0_decayTime/F");
184  Background_tree->Branch("outKFP_D0_decayLength", &outKFP_D0_decayLength,"outKFP_D0_decayLength/F");
185  Background_tree->Branch("outKFP_D0_DIRA", &outKFP_D0_DIRA,"outKFP_D0_DIRA/F");
186  Background_tree->Branch("outKFP_D0_IP", &outKFP_D0_IP,"outKFP_D0_IP/F");
187  Background_tree->Branch("outKFP_D0_IPchi2", &outKFP_D0_IPchi2,"outKFP_D0_IPchi2/F");
188  Background_tree->Branch("outKFP_D0_pseudorapidity", &outKFP_D0_pseudorapidity,"outKFP_D0_pseudorapidity/F");
189  Background_tree->Branch("outKFP_D0_rapidity", &outKFP_D0_rapidity,"outKFP_D0_rapidity/F");
190  Background_tree->Branch("outKFP_positive_px", &outKFP_positive_px,"outKFP_positive_px/F");
191  Background_tree->Branch("outKFP_positive_py", &outKFP_positive_py,"outKFP_positive_py/F");
192  Background_tree->Branch("outKFP_positive_pz", &outKFP_positive_pz,"outKFP_positive_pz/F");
193  Background_tree->Branch("outKFP_positive_pE", &outKFP_positive_pE,"outKFP_positive_pE/F");
194  Background_tree->Branch("outKFP_negative_px", &outKFP_negative_px,"outKFP_negative_px/F");
195  Background_tree->Branch("outKFP_negative_py", &outKFP_negative_py,"outKFP_negative_py/F");
196  Background_tree->Branch("outKFP_negative_pz", &outKFP_negative_pz,"outKFP_negative_pz/F");
197  Background_tree->Branch("outKFP_negative_pE", &outKFP_negative_pE,"outKFP_negative_pE/F");
198  Background_tree->Branch("outKFP_positive_p", &outKFP_positive_p,"outKFP_positive_p/F");
199  Background_tree->Branch("outKFP_negative_p", &outKFP_negative_p,"outKFP_negative_p/F");
200  Background_tree->Branch("outKFP_positive_pseudorapidity", &outKFP_positive_pseudorapidity,"outKFP_positive_pseudorapidity/F");
201  Background_tree->Branch("outKFP_negative_pseudorapidity", &outKFP_negative_pseudorapidity,"outKFP_negative_pseudorapidity/F");
202  Background_tree->Branch("outKFP_positive_charge", &outKFP_positive_charge,"outKFP_positive_charge/I");
203  Background_tree->Branch("outKFP_negative_charge", &outKFP_negative_charge,"outKFP_negative_charge/I");
204  Background_tree->Branch("outKFP_positive_IP", &outKFP_positive_IP,"outKFP_positive_IP/F");
205  Background_tree->Branch("outKFP_negative_IP", &outKFP_negative_IP,"outKFP_negative_IP/F");
206  Background_tree->Branch("outKFP_positive_IPchi2", &outKFP_positive_IPchi2,"outKFP_positive_IPchi2/F");
207  Background_tree->Branch("outKFP_negative_IPchi2", &outKFP_negative_IPchi2,"outKFP_negative_IPchi2/F");
208  Background_tree->Branch("outKFP_track_1_track_2_DCA", &outKFP_track_1_track_2_DCA,"outKFP_track_1_track_2_DCA/F");
209  Background_tree->Branch("outKFP_KpPm_invm", &outKFP_KpPm_invm,"outKFP_KpPm_invm/F");
210  Background_tree->Branch("outKFP_KmPp_invm", &outKFP_KmPp_invm,"outKFP_KmPp_invm/F");
211 
212  for (int i = 1; i < 237; ++i)
213  {
214  std::cout << "Starting File Number: " << i << std::endl;
215  string fNumb;
216  if (i < 10)
217  {
218  fNumb = "00" + std::to_string(i);
219  }
220  else if (i < 100)
221  {
222  fNumb = "0" + std::to_string(i);
223  }
224  else
225  {
226  fNumb = std::to_string(i);
227  }
228 
229  string truthFile = "outputData_Run40Acceptance082922_00" + fNumb + ".root";
230  string KFPFile = "outputData_Run40Acceptance082922_KFP_00" + fNumb + ".root";
231  truthFilePath = inDir + truthFile;
232  KFPFilePath = inDir + KFPFile;
233 
234  TFile *truthInput(0);
235  truthInput = TFile::Open(truthFilePath);
236  TTree *truthTree = (TTree*)truthInput->Get("QAG4SimulationTruthDecay");
237  TFile *KFPInput(0);
238  KFPInput = TFile::Open(KFPFilePath);
239  TTree *recoTree = (TTree*)KFPInput->Get("DecayTree");
240 
241  Float_t truth_mother_px, truth_mother_py, truth_mother_pz, truth_mother_pE;
242  Float_t truth_mother_pT, truth_mother_eta;
243  Int_t truth_mother_barcode, truth_mother_PDG_ID;
244  Float_t truth_track_1_px, truth_track_1_py, truth_track_1_pz, truth_track_1_pE;
245  Float_t truth_track_2_px, truth_track_2_py, truth_track_2_pz, truth_track_2_pE;
246  Float_t truth_track_1_eta, truth_track_2_eta;
247  Int_t truth_track_1_PDG_ID, truth_track_2_PDG_ID;
248  Float_t kfp_D0_px, kfp_D0_py, kfp_D0_pz, kfp_D0_pE;
249  Float_t kfp_D0_mass;
250  Float_t kfp_D0_decayTime, kfp_D0_decayLength, kfp_D0_DIRA, kfp_D0_IP, kfp_D0_IPchi2;
251  Float_t kfp_D0_pseudorapidity, kfp_D0_rapidity;
252  Float_t kfp_track_1_px, kfp_track_1_py, kfp_track_1_pz, kfp_track_1_pE;
253  Float_t kfp_track_2_px, kfp_track_2_py, kfp_track_2_pz, kfp_track_2_pE;
254  Float_t kfp_track_1_p, kfp_track_2_p;
255  Float_t kfp_track_1_true_px, kfp_track_1_true_py, kfp_track_1_true_pz;
256  Float_t kfp_track_2_true_px, kfp_track_2_true_py, kfp_track_2_true_pz;
257  Float_t kfp_track_1_pseudorapidity, kfp_track_2_pseudorapidity;
258  Int_t kfp_track_1_charge, kfp_track_2_charge;
259  Float_t kfp_track_1_IP, kfp_track_2_IP;
260  Float_t kfp_track_1_IPchi2, kfp_track_2_IPchi2;
261  Float_t kfp_track_1_track_2_DCA;
262 
263  truthTree->SetBranchAddress("mother_px", &truth_mother_px);
264  truthTree->SetBranchAddress("mother_py", &truth_mother_py);
265  truthTree->SetBranchAddress("mother_pz", &truth_mother_pz);
266  truthTree->SetBranchAddress("mother_pE", &truth_mother_pE);
267  truthTree->SetBranchAddress("mother_pT", &truth_mother_pT);
268  truthTree->SetBranchAddress("mother_eta", &truth_mother_eta);
269  truthTree->SetBranchAddress("mother_barcode", &truth_mother_barcode);
270  truthTree->SetBranchAddress("mother_PDG_ID", &truth_mother_PDG_ID);
271  truthTree->SetBranchAddress("track_1_px", &truth_track_1_px);
272  truthTree->SetBranchAddress("track_1_py", &truth_track_1_py);
273  truthTree->SetBranchAddress("track_1_pz", &truth_track_1_pz);
274  truthTree->SetBranchAddress("track_1_pE", &truth_track_1_pE);
275  truthTree->SetBranchAddress("track_2_px", &truth_track_2_px);
276  truthTree->SetBranchAddress("track_2_py", &truth_track_2_py);
277  truthTree->SetBranchAddress("track_2_pz", &truth_track_2_pz);
278  truthTree->SetBranchAddress("track_2_pE", &truth_track_2_pE);
279  truthTree->SetBranchAddress("track_1_eta", &truth_track_1_eta);
280  truthTree->SetBranchAddress("track_2_eta", &truth_track_2_eta);
281  truthTree->SetBranchAddress("track_1_PDG_ID", &truth_track_1_PDG_ID);
282  truthTree->SetBranchAddress("track_2_PDG_ID", &truth_track_2_PDG_ID);
283 
284  recoTree->SetBranchAddress("D0_px", &kfp_D0_px);
285  recoTree->SetBranchAddress("D0_py", &kfp_D0_py);
286  recoTree->SetBranchAddress("D0_pz", &kfp_D0_pz);
287  recoTree->SetBranchAddress("D0_pE", &kfp_D0_pE);
288  recoTree->SetBranchAddress("D0_mass", &kfp_D0_mass);
289  recoTree->SetBranchAddress("D0_decayTime", &kfp_D0_decayTime);
290  recoTree->SetBranchAddress("D0_decayLength", &kfp_D0_decayLength);
291  recoTree->SetBranchAddress("D0_DIRA", &kfp_D0_DIRA);
292  recoTree->SetBranchAddress("D0_IP", &kfp_D0_IP);
293  recoTree->SetBranchAddress("D0_IPchi2", &kfp_D0_IPchi2);
294  recoTree->SetBranchAddress("D0_pseudorapidity", &kfp_D0_pseudorapidity);
295  recoTree->SetBranchAddress("D0_rapidity", &kfp_D0_rapidity);
296  recoTree->SetBranchAddress("track_1_px", &kfp_track_1_px);
297  recoTree->SetBranchAddress("track_1_py", &kfp_track_1_py);
298  recoTree->SetBranchAddress("track_1_pz", &kfp_track_1_pz);
299  recoTree->SetBranchAddress("track_1_pE", &kfp_track_1_pE);
300  recoTree->SetBranchAddress("track_2_px", &kfp_track_2_px);
301  recoTree->SetBranchAddress("track_2_py", &kfp_track_2_py);
302  recoTree->SetBranchAddress("track_2_pz", &kfp_track_2_pz);
303  recoTree->SetBranchAddress("track_2_pE", &kfp_track_2_pE);
304  recoTree->SetBranchAddress("track_1_p", &kfp_track_1_p);
305  recoTree->SetBranchAddress("track_2_p", &kfp_track_2_p);
306  recoTree->SetBranchAddress("track_1_true_px", &kfp_track_1_true_px);
307  recoTree->SetBranchAddress("track_1_true_py", &kfp_track_1_true_py);
308  recoTree->SetBranchAddress("track_1_true_pz", &kfp_track_1_true_pz);
309  recoTree->SetBranchAddress("track_2_true_px", &kfp_track_2_true_px);
310  recoTree->SetBranchAddress("track_2_true_py", &kfp_track_2_true_py);
311  recoTree->SetBranchAddress("track_2_true_pz", &kfp_track_2_true_pz);
312  recoTree->SetBranchAddress("track_1_pseudorapidity", &kfp_track_1_pseudorapidity);
313  recoTree->SetBranchAddress("track_2_pseudorapidity", &kfp_track_2_pseudorapidity);
314  recoTree->SetBranchAddress("track_1_charge", &kfp_track_1_charge);
315  recoTree->SetBranchAddress("track_2_charge", &kfp_track_2_charge);
316  recoTree->SetBranchAddress("track_1_IP", &kfp_track_1_IP);
317  recoTree->SetBranchAddress("track_2_IP", &kfp_track_2_IP);
318  recoTree->SetBranchAddress("track_1_IPchi2", &kfp_track_1_IPchi2);
319  recoTree->SetBranchAddress("track_2_IPchi2", &kfp_track_2_IPchi2);
320  recoTree->SetBranchAddress("track_1_track_2_DCA", &kfp_track_1_track_2_DCA);
321 
322  vector<Int_t> reconstructedBarcodes;
323  vector<Int_t> missedBarcodes;
324  vector<int> usedRecoEntries;
325 
326  int minEntry;
327 
328  for (int j = 0; j < truthTree->GetEntries(); ++j)
329  {
330  truthTree->GetEntry(j);
331 
332  outTrue_mother_px = 0; outTrue_mother_py = 0; outTrue_mother_pz = 0; outTrue_mother_pE = 0;
333  outTrue_mother_pT = 0; outTrue_mother_eta = 0;
334  outTrue_mother_barcode = 0; outTrue_mother_PDG_ID = 0;
335  outTrue_positive_px = 0; outTrue_positive_py = 0; outTrue_positive_pz = 0; outTrue_positive_pE = 0;
336  outTrue_negative_px = 0; outTrue_negative_py = 0; outTrue_negative_pz = 0; outTrue_negative_pE = 0;
337  outTrue_positive_eta = 0; outTrue_negative_eta = 0;
338  outTrue_positive_PDG_ID = 0; outTrue_negative_PDG_ID = 0;
339  outKFP_D0_px = 0; outKFP_D0_py = 0; outKFP_D0_pz = 0; outKFP_D0_pE = 0;
340  outKFP_D0_mass = 0;
341  outKFP_D0_decayTime = 0; outKFP_D0_decayLength = 0; outKFP_D0_DIRA = 0; outKFP_D0_IP = 0; outKFP_D0_IPchi2 = 0;
342  outKFP_D0_pseudorapidity = 0; outKFP_D0_rapidity = 0;
343  outKFP_positive_px = 0; outKFP_positive_py = 0; outKFP_positive_pz = 0; outKFP_positive_pE = 0;
344  outKFP_negative_px = 0; outKFP_negative_py = 0; outKFP_negative_pz = 0; outKFP_negative_pE = 0;
345  outKFP_positive_p = 0; outKFP_negative_p = 0;
346  outKFP_positive_pseudorapidity = 0; outKFP_negative_pseudorapidity = 0;
347  outKFP_positive_charge = 0; outKFP_negative_charge = 0;
348  outKFP_positive_IP = 0; outKFP_negative_IP = 0;
349  outKFP_positive_IPchi2 = 0; outKFP_negative_IPchi2 = 0;
350  outKFP_track_1_track_2_DCA = 0;
351  outKFP_KpPm_invm = 0; outKFP_KmPp_invm = 0;
352 
353  bool matchPx_t1r1;
354  bool matchPy_t1r1;
355  bool matchPz_t1r1;
356  bool matchPx_t1r2;
357  bool matchPy_t1r2;
358  bool matchPz_t1r2;
359  bool matchPx_t2r1;
360  bool matchPy_t2r1;
361  bool matchPz_t2r1;
362  bool matchPx_t2r2;
363  bool matchPy_t2r2;
364  bool matchPz_t2r2;
365 
366  bool match_t1r1_t2r2 = false;
367  bool match_t2r1_t1r2 = false;
368 
369  for (int k = 0; k < recoTree->GetEntries(); ++k)
370  {
371  recoTree->GetEntry(k);
372 
373  matchPx_t1r1 = (kfp_track_1_true_px == truth_track_1_px);
374  matchPy_t1r1 = (kfp_track_1_true_py == truth_track_1_py);
375  matchPz_t1r1 = (kfp_track_1_true_pz == truth_track_1_pz);
376  matchPx_t1r2 = (kfp_track_2_true_px == truth_track_1_px);
377  matchPy_t1r2 = (kfp_track_2_true_py == truth_track_1_py);
378  matchPz_t1r2 = (kfp_track_2_true_pz == truth_track_1_pz);
379  matchPx_t2r1 = (kfp_track_1_true_px == truth_track_2_px);
380  matchPy_t2r1 = (kfp_track_1_true_py == truth_track_2_py);
381  matchPz_t2r1 = (kfp_track_1_true_pz == truth_track_2_pz);
382  matchPx_t2r2 = (kfp_track_2_true_px == truth_track_2_px);
383  matchPy_t2r2 = (kfp_track_2_true_py == truth_track_2_py);
384  matchPz_t2r2 = (kfp_track_2_true_pz == truth_track_2_pz);
385 
386  if (matchPx_t1r1 && matchPy_t1r1 && matchPz_t1r1)
387  {
388  if (matchPx_t2r2 && matchPy_t2r2 && matchPz_t2r2)
389  {
390  std::cout << "Got a match" << std::endl;
391  match_t1r1_t2r2 = true;
392  minEntry = k;
393  break;
394  }
395  }
396  else if (matchPx_t1r2 && matchPy_t1r2 && matchPz_t1r2)
397  {
398  if (matchPx_t2r1 && matchPy_t2r1 && matchPz_t2r1)
399  {
400  std::cout << "Got a match" << std::endl;
401  match_t2r1_t1r2 = true;
402  minEntry = k;
403  break;
404  }
405  }
406  }
407  if (match_t1r1_t2r2)
408  {
409  if (kfp_track_1_charge > 0 && kfp_track_2_charge < 0)
410  {
411  outTrue_mother_px = truth_mother_px; outTrue_mother_py = truth_mother_py; outTrue_mother_pz = truth_mother_pz; outTrue_mother_pE = truth_mother_pE;
412  outTrue_mother_pT = truth_mother_pT; outTrue_mother_eta = truth_mother_eta;
413  outTrue_mother_barcode = truth_mother_barcode; outTrue_mother_PDG_ID = truth_mother_PDG_ID;
414  outTrue_positive_px = truth_track_1_px; outTrue_positive_py = truth_track_1_py; outTrue_positive_pz = truth_track_1_pz; outTrue_positive_pE = truth_track_1_pE;
415  outTrue_negative_px = truth_track_2_px; outTrue_negative_py = truth_track_2_py; outTrue_negative_pz = truth_track_2_pz; outTrue_negative_pE = truth_track_2_pE;
416  outTrue_positive_eta = truth_track_1_eta; outTrue_negative_eta = truth_track_2_eta;
417  outTrue_positive_PDG_ID = truth_track_1_PDG_ID; outTrue_negative_PDG_ID = truth_track_2_PDG_ID;
418  outKFP_D0_px = kfp_D0_px; outKFP_D0_py = kfp_D0_py; outKFP_D0_pz = kfp_D0_pz; outKFP_D0_pE = kfp_D0_pE;
419  outKFP_D0_mass = kfp_D0_mass;
420  outKFP_D0_decayTime = kfp_D0_decayTime; outKFP_D0_decayLength = kfp_D0_decayLength; outKFP_D0_DIRA = kfp_D0_DIRA; outKFP_D0_IP = kfp_D0_IP; outKFP_D0_IPchi2 = kfp_D0_IPchi2;
421  outKFP_D0_pseudorapidity = kfp_D0_pseudorapidity; outKFP_D0_rapidity = kfp_D0_rapidity;
422  outKFP_positive_px = kfp_track_1_px; outKFP_positive_py = kfp_track_1_py; outKFP_positive_pz = kfp_track_1_pz; outKFP_positive_pE = kfp_track_1_pE;
423  outKFP_negative_px = kfp_track_2_px; outKFP_negative_py = kfp_track_2_py; outKFP_negative_pz = kfp_track_2_pz; outKFP_negative_pE = kfp_track_2_pE;
424  outKFP_positive_p = kfp_track_1_p; outKFP_negative_p = kfp_track_2_p;
425  outKFP_positive_pseudorapidity = kfp_track_1_pseudorapidity; outKFP_negative_pseudorapidity = kfp_track_2_pseudorapidity;
426  outKFP_positive_charge = kfp_track_1_charge; outKFP_negative_charge = kfp_track_2_charge;
427  outKFP_positive_IP = kfp_track_1_IP; outKFP_negative_IP = kfp_track_2_IP;
428  outKFP_positive_IPchi2 = kfp_track_1_IPchi2; outKFP_negative_IPchi2 = kfp_track_2_IPchi2;
429  outKFP_track_1_track_2_DCA = kfp_track_1_track_2_DCA;
430 
431  TVector3 P_p(outKFP_positive_px,outKFP_positive_py,outKFP_positive_pz);
432  TVector3 N_p(outKFP_negative_px,outKFP_negative_py,outKFP_negative_pz);
433  Float_t K_energy = sqrt(P_p.Mag2() + pow(K_mass,2));
434  Float_t Pi_energy = sqrt(N_p.Mag2() + pow(Pi_mass,2));
435  outKFP_KpPm_invm = sqrt(pow((K_energy + Pi_energy),2) - ((N_p+P_p).Mag2()));
436  K_energy = sqrt(N_p.Mag2() + pow(K_mass,2));
437  Pi_energy = sqrt(P_p.Mag2() + pow(Pi_mass,2));
438  outKFP_KmPp_invm = sqrt(pow((K_energy + Pi_energy),2) - ((N_p+P_p).Mag2()));
439 
440  if (usedRecoEntries.size() > 0)
441  {
442  bool recoCandidateUsed = false;
443  for (int ent : usedRecoEntries)
444  {
445  if (ent == minEntry)
446  {
447  std::cout << "Candidate has already been used, skipping this candidate" << std::endl;
448  recoCandidateUsed = true;
449  }
450  }
451  if (recoCandidateUsed) continue;
452  else
453  {
454  usedRecoEntries.push_back(minEntry);
455  // D0 or D0bar check here
456  if (outTrue_mother_PDG_ID == 421) D0_tree->Fill();
457  else if (outTrue_mother_PDG_ID == -421) D0bar_tree->Fill();
458  else
459  {
460  Background_tree->Fill();
461  std::cout << "ERROR: Should not get here if tracks match" <<std::endl;
462  }
463  }
464  }
465  else
466  {
467  usedRecoEntries.push_back(minEntry);
468 
469  // D0 or D0bar check here
470  if (outTrue_mother_PDG_ID == 421) D0_tree->Fill();
471  else if (outTrue_mother_PDG_ID == -421) D0bar_tree->Fill();
472  else
473  {
474  Background_tree->Fill();
475  std::cout << "ERROR: Should not get here if tracks match" <<std::endl;
476  }
477  }
478  }
479  else if (kfp_track_1_charge < 0 && kfp_track_2_charge > 0)
480  {
481  outTrue_mother_px = truth_mother_px; outTrue_mother_py = truth_mother_py; outTrue_mother_pz = truth_mother_pz; outTrue_mother_pE = truth_mother_pE;
482  outTrue_mother_pT = truth_mother_pT; outTrue_mother_eta = truth_mother_eta;
483  outTrue_mother_barcode = truth_mother_barcode; outTrue_mother_PDG_ID = truth_mother_PDG_ID;
484  outTrue_negative_px = truth_track_1_px; outTrue_negative_py = truth_track_1_py; outTrue_negative_pz = truth_track_1_pz; outTrue_negative_pE = truth_track_1_pE;
485  outTrue_positive_px = truth_track_2_px; outTrue_positive_py = truth_track_2_py; outTrue_positive_pz = truth_track_2_pz; outTrue_positive_pE = truth_track_2_pE;
486  outTrue_negative_eta = truth_track_1_eta; outTrue_positive_eta = truth_track_2_eta;
487  outTrue_negative_PDG_ID = truth_track_1_PDG_ID; outTrue_positive_PDG_ID = truth_track_2_PDG_ID;
488  outKFP_D0_px = kfp_D0_px; outKFP_D0_py = kfp_D0_py; outKFP_D0_pz = kfp_D0_pz; outKFP_D0_pE = kfp_D0_pE;
489  outKFP_D0_mass = kfp_D0_mass;
490  outKFP_D0_decayTime = kfp_D0_decayTime; outKFP_D0_decayLength = kfp_D0_decayLength; outKFP_D0_DIRA = kfp_D0_DIRA; outKFP_D0_IP = kfp_D0_IP; outKFP_D0_IPchi2 = kfp_D0_IPchi2;
491  outKFP_D0_pseudorapidity = kfp_D0_pseudorapidity; outKFP_D0_rapidity = kfp_D0_rapidity;
492  outKFP_negative_px = kfp_track_1_px; outKFP_negative_py = kfp_track_1_py; outKFP_negative_pz = kfp_track_1_pz; outKFP_negative_pE = kfp_track_1_pE;
493  outKFP_positive_px = kfp_track_2_px; outKFP_positive_py = kfp_track_2_py; outKFP_positive_pz = kfp_track_2_pz; outKFP_positive_pE = kfp_track_2_pE;
494  outKFP_negative_p = kfp_track_1_p; outKFP_positive_p = kfp_track_2_p;
495  outKFP_negative_pseudorapidity = kfp_track_1_pseudorapidity; outKFP_positive_pseudorapidity = kfp_track_2_pseudorapidity;
496  outKFP_negative_charge = kfp_track_1_charge; outKFP_positive_charge = kfp_track_2_charge;
497  outKFP_negative_IP = kfp_track_1_IP; outKFP_positive_IP = kfp_track_2_IP;
498  outKFP_negative_IPchi2 = kfp_track_1_IPchi2; outKFP_positive_IPchi2 = kfp_track_2_IPchi2;
499  outKFP_track_1_track_2_DCA = kfp_track_1_track_2_DCA;
500 
501  TVector3 P_p(outKFP_positive_px,outKFP_positive_py,outKFP_positive_pz);
502  TVector3 N_p(outKFP_negative_px,outKFP_negative_py,outKFP_negative_pz);
503  Float_t K_energy = sqrt(P_p.Mag2() + pow(K_mass,2));
504  Float_t Pi_energy = sqrt(N_p.Mag2() + pow(Pi_mass,2));
505  outKFP_KpPm_invm = sqrt(pow((K_energy + Pi_energy),2) - ((N_p+P_p).Mag2()));
506  K_energy = sqrt(N_p.Mag2() + pow(K_mass,2));
507  Pi_energy = sqrt(P_p.Mag2() + pow(Pi_mass,2));
508  outKFP_KmPp_invm = sqrt(pow((K_energy + Pi_energy),2) - ((N_p+P_p).Mag2()));
509 
510  if (usedRecoEntries.size() > 0)
511  {
512  bool recoCandidateUsed = false;
513  for (int ent : usedRecoEntries)
514  {
515  if (ent == minEntry)
516  {
517  std::cout << "Candidate has already been used, skipping this candidate" << std::endl;
518  recoCandidateUsed = true;
519  }
520  }
521  if (recoCandidateUsed) continue;
522  else
523  {
524  usedRecoEntries.push_back(minEntry);
525  // D0 or D0bar check here
526  if (outTrue_mother_PDG_ID == 421) D0_tree->Fill();
527  else if (outTrue_mother_PDG_ID == -421) D0bar_tree->Fill();
528  else
529  {
530  Background_tree->Fill();
531  std::cout << "ERROR: Should not get here if tracks match" <<std::endl;
532  }
533  }
534  }
535  else
536  {
537  usedRecoEntries.push_back(minEntry);
538 
539  // D0 or D0bar check here
540  if (outTrue_mother_PDG_ID == 421) D0_tree->Fill();
541  else if (outTrue_mother_PDG_ID == -421) D0bar_tree->Fill();
542  else
543  {
544  Background_tree->Fill();
545  std::cout << "ERROR: Should not get here if tracks match" <<std::endl;
546  }
547  }
548  }
549  }
550  else if (match_t2r1_t1r2)
551  {
552  if (kfp_track_1_charge > 0 && kfp_track_2_charge < 0)
553  {
554  outTrue_mother_px = truth_mother_px; outTrue_mother_py = truth_mother_py; outTrue_mother_pz = truth_mother_pz; outTrue_mother_pE = truth_mother_pE;
555  outTrue_mother_pT = truth_mother_pT; outTrue_mother_eta = truth_mother_eta;
556  outTrue_mother_barcode = truth_mother_barcode; outTrue_mother_PDG_ID = truth_mother_PDG_ID;
557  outTrue_negative_px = truth_track_1_px; outTrue_negative_py = truth_track_1_py; outTrue_negative_pz = truth_track_1_pz; outTrue_negative_pE = truth_track_1_pE;
558  outTrue_positive_px = truth_track_2_px; outTrue_positive_py = truth_track_2_py; outTrue_positive_pz = truth_track_2_pz; outTrue_positive_pE = truth_track_2_pE;
559  outTrue_negative_eta = truth_track_1_eta; outTrue_positive_eta = truth_track_2_eta;
560  outTrue_negative_PDG_ID = truth_track_1_PDG_ID; outTrue_positive_PDG_ID = truth_track_2_PDG_ID;
561  outKFP_D0_px = kfp_D0_px; outKFP_D0_py = kfp_D0_py; outKFP_D0_pz = kfp_D0_pz; outKFP_D0_pE = kfp_D0_pE;
562  outKFP_D0_mass = kfp_D0_mass;
563  outKFP_D0_decayTime = kfp_D0_decayTime; outKFP_D0_decayLength = kfp_D0_decayLength; outKFP_D0_DIRA = kfp_D0_DIRA; outKFP_D0_IP = kfp_D0_IP; outKFP_D0_IPchi2 = kfp_D0_IPchi2;
564  outKFP_D0_pseudorapidity = kfp_D0_pseudorapidity; outKFP_D0_rapidity = kfp_D0_rapidity;
565  outKFP_positive_px = kfp_track_1_px; outKFP_positive_py = kfp_track_1_py; outKFP_positive_pz = kfp_track_1_pz; outKFP_positive_pE = kfp_track_1_pE;
566  outKFP_negative_px = kfp_track_2_px; outKFP_negative_py = kfp_track_2_py; outKFP_negative_pz = kfp_track_2_pz; outKFP_negative_pE = kfp_track_2_pE;
567  outKFP_positive_p = kfp_track_1_p; outKFP_negative_p = kfp_track_2_p;
568  outKFP_positive_pseudorapidity = kfp_track_1_pseudorapidity; outKFP_negative_pseudorapidity = kfp_track_2_pseudorapidity;
569  outKFP_positive_charge = kfp_track_1_charge; outKFP_negative_charge = kfp_track_2_charge;
570  outKFP_positive_IP = kfp_track_1_IP; outKFP_negative_IP = kfp_track_2_IP;
571  outKFP_positive_IPchi2 = kfp_track_1_IPchi2; outKFP_negative_IPchi2 = kfp_track_2_IPchi2;
572  outKFP_track_1_track_2_DCA = kfp_track_1_track_2_DCA;
573 
574  TVector3 P_p(outKFP_positive_px,outKFP_positive_py,outKFP_positive_pz);
575  TVector3 N_p(outKFP_negative_px,outKFP_negative_py,outKFP_negative_pz);
576  Float_t K_energy = sqrt(P_p.Mag2() + pow(K_mass,2));
577  Float_t Pi_energy = sqrt(N_p.Mag2() + pow(Pi_mass,2));
578  outKFP_KpPm_invm = sqrt(pow((K_energy + Pi_energy),2) - ((N_p+P_p).Mag2()));
579  K_energy = sqrt(N_p.Mag2() + pow(K_mass,2));
580  Pi_energy = sqrt(P_p.Mag2() + pow(Pi_mass,2));
581  outKFP_KmPp_invm = sqrt(pow((K_energy + Pi_energy),2) - ((N_p+P_p).Mag2()));
582 
583  if (usedRecoEntries.size() > 0)
584  {
585  bool recoCandidateUsed = false;
586  for (int ent : usedRecoEntries)
587  {
588  if (ent == minEntry)
589  {
590  std::cout << "Candidate has already been used, skipping this candidate" << std::endl;
591  recoCandidateUsed = true;
592  }
593  }
594  if (recoCandidateUsed) continue;
595  else
596  {
597  usedRecoEntries.push_back(minEntry);
598  // D0 or D0bar check here
599  if (outTrue_mother_PDG_ID == 421) D0_tree->Fill();
600  else if (outTrue_mother_PDG_ID == -421) D0bar_tree->Fill();
601  else
602  {
603  Background_tree->Fill();
604  std::cout << "ERROR: Should not get here if tracks match" <<std::endl;
605  }
606  }
607  }
608  else
609  {
610  usedRecoEntries.push_back(minEntry);
611 
612  // D0 or D0bar check here
613  if (outTrue_mother_PDG_ID == 421) D0_tree->Fill();
614  else if (outTrue_mother_PDG_ID == -421) D0bar_tree->Fill();
615  else
616  {
617  Background_tree->Fill();
618  std::cout << "ERROR: Should not get here if tracks match" <<std::endl;
619  }
620  }
621  }
622  if (kfp_track_1_charge < 0 && kfp_track_2_charge > 0)
623  {
624  outTrue_mother_px = truth_mother_px; outTrue_mother_py = truth_mother_py; outTrue_mother_pz = truth_mother_pz; outTrue_mother_pE = truth_mother_pE;
625  outTrue_mother_pT = truth_mother_pT; outTrue_mother_eta = truth_mother_eta;
626  outTrue_mother_barcode = truth_mother_barcode; outTrue_mother_PDG_ID = truth_mother_PDG_ID;
627  outTrue_positive_px = truth_track_1_px; outTrue_positive_py = truth_track_1_py; outTrue_positive_pz = truth_track_1_pz; outTrue_positive_pE = truth_track_1_pE;
628  outTrue_negative_px = truth_track_2_px; outTrue_negative_py = truth_track_2_py; outTrue_negative_pz = truth_track_2_pz; outTrue_negative_pE = truth_track_2_pE;
629  outTrue_positive_eta = truth_track_1_eta; outTrue_negative_eta = truth_track_2_eta;
630  outTrue_positive_PDG_ID = truth_track_1_PDG_ID; outTrue_negative_PDG_ID = truth_track_2_PDG_ID;
631  outKFP_D0_px = kfp_D0_px; outKFP_D0_py = kfp_D0_py; outKFP_D0_pz = kfp_D0_pz; outKFP_D0_pE = kfp_D0_pE;
632  outKFP_D0_mass = kfp_D0_mass;
633  outKFP_D0_decayTime = kfp_D0_decayTime; outKFP_D0_decayLength = kfp_D0_decayLength; outKFP_D0_DIRA = kfp_D0_DIRA; outKFP_D0_IP = kfp_D0_IP; outKFP_D0_IPchi2 = kfp_D0_IPchi2;
634  outKFP_D0_pseudorapidity = kfp_D0_pseudorapidity; outKFP_D0_rapidity = kfp_D0_rapidity;
635  outKFP_negative_px = kfp_track_1_px; outKFP_negative_py = kfp_track_1_py; outKFP_negative_pz = kfp_track_1_pz; outKFP_negative_pE = kfp_track_1_pE;
636  outKFP_positive_px = kfp_track_2_px; outKFP_positive_py = kfp_track_2_py; outKFP_positive_pz = kfp_track_2_pz; outKFP_positive_pE = kfp_track_2_pE;
637  outKFP_negative_p = kfp_track_1_p; outKFP_positive_p = kfp_track_2_p;
638  outKFP_negative_pseudorapidity = kfp_track_1_pseudorapidity; outKFP_positive_pseudorapidity = kfp_track_2_pseudorapidity;
639  outKFP_negative_charge = kfp_track_1_charge; outKFP_positive_charge = kfp_track_2_charge;
640  outKFP_negative_IP = kfp_track_1_IP; outKFP_positive_IP = kfp_track_2_IP;
641  outKFP_negative_IPchi2 = kfp_track_1_IPchi2; outKFP_positive_IPchi2 = kfp_track_2_IPchi2;
642  outKFP_track_1_track_2_DCA = kfp_track_1_track_2_DCA;
643 
644  TVector3 P_p(outKFP_positive_px,outKFP_positive_py,outKFP_positive_pz);
645  TVector3 N_p(outKFP_negative_px,outKFP_negative_py,outKFP_negative_pz);
646  Float_t K_energy = sqrt(P_p.Mag2() + pow(K_mass,2));
647  Float_t Pi_energy = sqrt(N_p.Mag2() + pow(Pi_mass,2));
648  outKFP_KpPm_invm = sqrt(pow((K_energy + Pi_energy),2) - ((N_p+P_p).Mag2()));
649  K_energy = sqrt(N_p.Mag2() + pow(K_mass,2));
650  Pi_energy = sqrt(P_p.Mag2() + pow(Pi_mass,2));
651  outKFP_KmPp_invm = sqrt(pow((K_energy + Pi_energy),2) - ((N_p+P_p).Mag2()));
652 
653  if (usedRecoEntries.size() > 0)
654  {
655  bool recoCandidateUsed = false;
656  for (int ent : usedRecoEntries)
657  {
658  if (ent == minEntry)
659  {
660  std::cout << "Candidate has already been used, skipping this candidate" << std::endl;
661  recoCandidateUsed = true;
662  }
663  }
664  if (recoCandidateUsed) continue;
665  else
666  {
667  usedRecoEntries.push_back(minEntry);
668  // D0 or D0bar check here
669  if (outTrue_mother_PDG_ID == 421) D0_tree->Fill();
670  else if (outTrue_mother_PDG_ID == -421) D0bar_tree->Fill();
671  else
672  {
673  Background_tree->Fill();
674  std::cout << "ERROR: Should not get here if tracks match" <<std::endl;
675  }
676  }
677  }
678  else
679  {
680  usedRecoEntries.push_back(minEntry);
681 
682  // D0 or D0bar check here
683  if (outTrue_mother_PDG_ID == 421) D0_tree->Fill();
684  else if (outTrue_mother_PDG_ID == -421) D0bar_tree->Fill();
685  else
686  {
687  Background_tree->Fill();
688  std::cout << "ERROR: Should not get here if tracks match" <<std::endl;
689  }
690  }
691  }
692  }
693  }
694  for (int k = 0; k < recoTree->GetEntries(); ++k)
695  {
696  bool recoCandidateUsed = false;
697  for (int ent : usedRecoEntries)
698  {
699  if (ent == k) recoCandidateUsed = true;
700  }
701  if (recoCandidateUsed == false)
702  {
703  recoTree->GetEntry(k);
704 
705  if (kfp_track_1_charge > 0 && kfp_track_2_charge < 0)
706  {
707  outKFP_D0_px = kfp_D0_px; outKFP_D0_py = kfp_D0_py; outKFP_D0_pz = kfp_D0_pz; outKFP_D0_pE = kfp_D0_pE;
708  outKFP_D0_mass = kfp_D0_mass;
709  outKFP_D0_decayTime = kfp_D0_decayTime; outKFP_D0_decayLength = kfp_D0_decayLength; outKFP_D0_DIRA = kfp_D0_DIRA; outKFP_D0_IP = kfp_D0_IP; outKFP_D0_IPchi2 = kfp_D0_IPchi2;
710  outKFP_D0_pseudorapidity = kfp_D0_pseudorapidity; outKFP_D0_rapidity = kfp_D0_rapidity;
711  outKFP_positive_px = kfp_track_1_px; outKFP_positive_py = kfp_track_1_py; outKFP_positive_pz = kfp_track_1_pz; outKFP_positive_pE = kfp_track_1_pE;
712  outKFP_negative_px = kfp_track_2_px; outKFP_negative_py = kfp_track_2_py; outKFP_negative_pz = kfp_track_2_pz; outKFP_negative_pE = kfp_track_2_pE;
713  outKFP_positive_p = kfp_track_1_p; outKFP_negative_p = kfp_track_2_p;
714  outKFP_positive_pseudorapidity = kfp_track_1_pseudorapidity; outKFP_negative_pseudorapidity = kfp_track_2_pseudorapidity;
715  outKFP_positive_charge = kfp_track_1_charge; outKFP_negative_charge = kfp_track_2_charge;
716  outKFP_positive_IP = kfp_track_1_IP; outKFP_negative_IP = kfp_track_2_IP;
717  outKFP_positive_IPchi2 = kfp_track_1_IPchi2; outKFP_negative_IPchi2 = kfp_track_2_IPchi2;
718  outKFP_track_1_track_2_DCA = kfp_track_1_track_2_DCA;
719 
720  TVector3 P_p(outKFP_positive_px,outKFP_positive_py,outKFP_positive_pz);
721  TVector3 N_p(outKFP_negative_px,outKFP_negative_py,outKFP_negative_pz);
722  Float_t K_energy = sqrt(P_p.Mag2() + pow(K_mass,2));
723  Float_t Pi_energy = sqrt(N_p.Mag2() + pow(Pi_mass,2));
724  outKFP_KpPm_invm = sqrt(pow((K_energy + Pi_energy),2) - ((N_p+P_p).Mag2()));
725  K_energy = sqrt(N_p.Mag2() + pow(K_mass,2));
726  Pi_energy = sqrt(P_p.Mag2() + pow(Pi_mass,2));
727  outKFP_KmPp_invm = sqrt(pow((K_energy + Pi_energy),2) - ((N_p+P_p).Mag2()));
728 
729  Background_tree->Fill();
730  }
731  else if (kfp_track_1_charge < 0 && kfp_track_2_charge > 0)
732  {
733  outKFP_D0_px = kfp_D0_px; outKFP_D0_py = kfp_D0_py; outKFP_D0_pz = kfp_D0_pz; outKFP_D0_pE = kfp_D0_pE;
734  outKFP_D0_mass = kfp_D0_mass;
735  outKFP_D0_decayTime = kfp_D0_decayTime; outKFP_D0_decayLength = kfp_D0_decayLength; outKFP_D0_DIRA = kfp_D0_DIRA; outKFP_D0_IP = kfp_D0_IP; outKFP_D0_IPchi2 = kfp_D0_IPchi2;
736  outKFP_D0_pseudorapidity = kfp_D0_pseudorapidity; outKFP_D0_rapidity = kfp_D0_rapidity;
737  outKFP_negative_px = kfp_track_1_px; outKFP_negative_py = kfp_track_1_py; outKFP_negative_pz = kfp_track_1_pz; outKFP_negative_pE = kfp_track_1_pE;
738  outKFP_positive_px = kfp_track_2_px; outKFP_positive_py = kfp_track_2_py; outKFP_positive_pz = kfp_track_2_pz; outKFP_positive_pE = kfp_track_2_pE;
739  outKFP_negative_p = kfp_track_1_p; outKFP_positive_p = kfp_track_2_p;
740  outKFP_negative_pseudorapidity = kfp_track_1_pseudorapidity; outKFP_positive_pseudorapidity = kfp_track_2_pseudorapidity;
741  outKFP_negative_charge = kfp_track_1_charge; outKFP_positive_charge = kfp_track_2_charge;
742  outKFP_negative_IP = kfp_track_1_IP; outKFP_positive_IP = kfp_track_2_IP;
743  outKFP_negative_IPchi2 = kfp_track_1_IPchi2; outKFP_positive_IPchi2 = kfp_track_2_IPchi2;
744  outKFP_track_1_track_2_DCA = kfp_track_1_track_2_DCA;
745 
746  TVector3 P_p(outKFP_positive_px,outKFP_positive_py,outKFP_positive_pz);
747  TVector3 N_p(outKFP_negative_px,outKFP_negative_py,outKFP_negative_pz);
748  Float_t K_energy = sqrt(P_p.Mag2() + pow(K_mass,2));
749  Float_t Pi_energy = sqrt(N_p.Mag2() + pow(Pi_mass,2));
750  outKFP_KpPm_invm = sqrt(pow((K_energy + Pi_energy),2) - ((N_p+P_p).Mag2()));
751  K_energy = sqrt(N_p.Mag2() + pow(K_mass,2));
752  Pi_energy = sqrt(P_p.Mag2() + pow(Pi_mass,2));
753  outKFP_KmPp_invm = sqrt(pow((K_energy + Pi_energy),2) - ((N_p+P_p).Mag2()));
754 
755  Background_tree->Fill();
756  }
757  else
758  {
759  std::cout << "WARNING: tracks not oppositely charged" << std::endl;
760  }
761  }
762  }
763  }
764 
765  newfile->cd();
766  D0_tree->Print();
767  D0_tree->Write();
768  D0bar_tree->Print();
769  D0bar_tree->Write();
770  Background_tree->Print();
771  Background_tree->Write();
772  newfile->Close();
773 
774  ifstream file2("Run40_D0_Separated_091922.root");
775  if (file2.good())
776  {
777  string moveOutput = "mv Run40_D0_Separated_091922.root " + inDir;
778  system(moveOutput.c_str());
779  }
780 
781 }