Changeset d77b51d in git for modules/IdentificationMap.cc
- Timestamp:
- Sep 29, 2015, 2:08:10 PM (9 years ago)
- Branches:
- ImprovedOutputFile, Timing, dual_readout, llp, master
- Children:
- a98c7ef
- Parents:
- d870fc5 (diff), 06ec139 (diff)
Note: this is a merge changeset, the changes displayed below correspond to the merge itself.
Use the(diff)
links above to see all the changes relative to each parent. - File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
modules/IdentificationMap.cc
rd870fc5 rd77b51d 69 69 void IdentificationMap::Init() 70 70 { 71 // read efficiency formula72 73 74 71 TMisIDMap::iterator itEfficiencyMap; 75 72 ExRootConfParam param; … … 87 84 formula->Compile(param[i*3 + 2].GetString()); 88 85 pdg = param[i*3].GetInt(); 89 fEfficiencyMap.insert(make_pair(pdg,make_pair(param[i*3 + 1].GetInt(),formula))); 90 91 // cout<<param[i*3].GetInt()<<","<<param[i*3+1].GetInt()<<","<<param[i*3 + 2].GetString()<<endl; 92 86 fEfficiencyMap.insert(make_pair(pdg, make_pair(param[i*3 + 1].GetInt(), formula))); 93 87 } 94 88 … … 100 94 formula->Compile("1.0"); 101 95 102 fEfficiencyMap.insert(make_pair(0, make_pair(0,formula)));96 fEfficiencyMap.insert(make_pair(0, make_pair(0, formula))); 103 97 } 104 98 … … 126 120 if(formula) delete formula; 127 121 } 128 129 122 } 130 123 … … 134 127 { 135 128 Candidate *candidate; 136 Double_t pt, eta, phi ;129 Double_t pt, eta, phi, e; 137 130 TMisIDMap::iterator itEfficiencyMap; 138 131 pair <TMisIDMap::iterator, TMisIDMap::iterator> range; 139 132 DelphesFormula *formula; 140 Int_t pdg In, pdgOut, charge;133 Int_t pdgCodeIn, pdgCodeOut, charge; 141 134 142 Double_t P, Pi; 143 144 // cout<<"------------ New Event ------------"<<endl; 135 Double_t p, r, total; 145 136 146 137 fItInputArray->Reset(); … … 152 143 phi = candidatePosition.Phi(); 153 144 pt = candidateMomentum.Pt(); 154 pdgIn = candidate->PID; 145 e = candidateMomentum.E(); 146 147 pdgCodeIn = candidate->PID; 155 148 charge = candidate->Charge; 156 149 157 // cout<<"------------ New Candidate ------------"<<endl;158 // cout<<candidate->PID<<" "<<pt<<","<<eta<<","<<phi<<endl;150 // first check that PID of this particle is specified in the map 151 // otherwise, look for PID = 0 159 152 160 P = 1.0;153 itEfficiencyMap = fEfficiencyMap.find(pdgCodeIn); 161 154 162 //first check that PID of this particle is specified in cfg, if not set look for PID=0 163 164 itEfficiencyMap = fEfficiencyMap.find(pdgIn); 165 166 range = fEfficiencyMap.equal_range(pdgIn); 167 if(range.first == range.second) range = fEfficiencyMap.equal_range(-pdgIn); 155 range = fEfficiencyMap.equal_range(pdgCodeIn); 156 if(range.first == range.second) range = fEfficiencyMap.equal_range(-pdgCodeIn); 168 157 if(range.first == range.second) range = fEfficiencyMap.equal_range(0); 169 158 170 //loop over submap for this pid 171 for (TMisIDMap::iterator it=range.first; it!=range.second; ++it) 159 r = gRandom->Uniform(); 160 total = 0.0; 161 162 // loop over sub-map for this PID 163 for(TMisIDMap::iterator it = range.first; it != range.second; ++it) 172 164 { 165 formula = (it->second).second; 166 pdgCodeOut = (it->second).first; 173 167 174 formula = (it->second).second; 175 pdgOut = (it->second).first; 168 p = formula->Eval(pt, eta, phi, e); 176 169 177 Pi = formula->Eval(pt, eta); 178 179 // check that sum of probabilities does not exceed 1. 180 P = (P - Pi)/P; 181 182 if( P < 0.0 ) continue; 183 else 170 if(total <= r && r < total + p) 184 171 { 185 186 //randomly assign a PID to particle according to map 187 Double_t rndm = gRandom->Uniform(); 188 189 if(rndm > P) 190 { 191 //change PID of particle 192 if(pdgOut != 0) candidate->PID = charge*pdgOut; 193 fOutputArray->Add(candidate); 194 break; 195 } 172 // change PID of particle 173 if(pdgCodeOut != 0) candidate->PID = charge*pdgCodeOut; 174 fOutputArray->Add(candidate); 175 break; 196 176 } 197 177 178 total += p; 198 179 } 199 200 } 201 180 } 202 181 } 203 182
Note:
See TracChangeset
for help on using the changeset viewer.