Changeset 21eab4f in git for modules/PileUpMerger.cc
- Timestamp:
- Aug 26, 2016, 5:17:14 PM (8 years ago)
- Parents:
- ec5e04b (diff), 94cacb6 (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. - git-author:
- AlexandreMertens <alexandre.mertens@…> (08/26/16 17:17:14)
- git-committer:
- GitHub <noreply@…> (08/26/16 17:17:14)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
modules/PileUpMerger.cc
rec5e04b r21eab4f 115 115 TDatabasePDG *pdg = TDatabasePDG::Instance(); 116 116 TParticlePDG *pdgParticle; 117 Int_t pid ;117 Int_t pid, nch, nvtx = -1; 118 118 Float_t x, y, z, t, vx, vy; 119 Float_t px, py, pz, e ;120 Double_t dz, dphi, dt ;119 Float_t px, py, pz, e, pt; 120 Double_t dz, dphi, dt, sumpt2, dz0, dt0; 121 121 Int_t numberOfEvents, event, numberOfParticles; 122 122 Long64_t allEntries, entry; … … 132 132 fFunction->GetRandom2(dz, dt); 133 133 134 dz0 = -1.0e6; 135 dt0 = -1.0e6; 136 134 137 dt *= c_light*1.0E3; // necessary in order to make t in mm/c 135 138 dz *= 1.0E3; // necessary in order to make z in mm 139 140 //cout<<dz<<","<<dt<<endl; 141 136 142 vx = 0.0; 137 143 vy = 0.0; 144 138 145 numberOfParticles = fInputArray->GetEntriesFast(); 146 nch = 0; 147 sumpt2 = 0.0; 148 149 factory = GetFactory(); 150 vertex = factory->NewCandidate(); 151 139 152 while((candidate = static_cast<Candidate*>(fItInputArray->Next()))) 140 153 { … … 143 156 z = candidate->Position.Z(); 144 157 t = candidate->Position.T(); 145 candidate->Position.SetZ(z + dz); 146 candidate->Position.SetT(t + dt); 158 pt = candidate->Momentum.Pt(); 159 160 // take postion and time from first stable particle 161 if (dz0 < -999999.0) 162 dz0 = z; 163 if (dt0 < -999999.0) 164 dt0 = t; 165 166 // cancel any possible offset in position and time the input file 167 candidate->Position.SetZ(z - dz0 + dz); 168 candidate->Position.SetT(t - dt0 + dt); 169 170 candidate->IsPU = 0; 171 147 172 fParticleOutputArray->Add(candidate); 173 174 if(TMath::Abs(candidate->Charge) > 1.0E-9) 175 { 176 nch++; 177 sumpt2 += pt*pt; 178 vertex->AddCandidate(candidate); 179 } 148 180 } 149 181 150 182 if(numberOfParticles > 0) 151 183 { 152 vx /= numberOfParticles;153 vy /= numberOfParticles;184 vx /= sumpt2; 185 vy /= sumpt2; 154 186 } 155 187 156 factory = GetFactory(); 157 158 vertex = factory->NewCandidate(); 188 nvtx++; 159 189 vertex->Position.SetXYZT(vx, vy, dz, dt); 190 vertex->ClusterIndex = nvtx; 191 vertex->ClusterNDF = nch; 192 vertex->SumPT2 = sumpt2; 193 vertex->GenSumPT2 = sumpt2; 160 194 fVertexOutputArray->Add(vertex); 161 195 … … 170 204 numberOfEvents = gRandom->Integer(2*fMeanPileUp + 1); 171 205 break; 206 case 2: 207 numberOfEvents = fMeanPileUp; 208 break; 172 209 default: 173 210 numberOfEvents = gRandom->Poisson(fMeanPileUp); … … 176 213 177 214 allEntries = fReader->GetEntries(); 215 178 216 179 217 for(event = 0; event < numberOfEvents; ++event) … … 198 236 vx = 0.0; 199 237 vy = 0.0; 238 200 239 numberOfParticles = 0; 240 sumpt2 = 0.0; 241 242 //factory = GetFactory(); 243 vertex = factory->NewCandidate(); 244 201 245 while(fReader->ReadParticle(pid, x, y, z, t, px, py, pz, e)) 202 246 { … … 215 259 candidate->Momentum.SetPxPyPzE(px, py, pz, e); 216 260 candidate->Momentum.RotateZ(dphi); 261 pt = candidate->Momentum.Pt(); 217 262 218 263 x -= fInputBeamSpotX; … … 224 269 vx += candidate->Position.X(); 225 270 vy += candidate->Position.Y(); 271 226 272 ++numberOfParticles; 273 if(TMath::Abs(candidate->Charge) > 1.0E-9) 274 { 275 nch++; 276 sumpt2 += pt*pt; 277 vertex->AddCandidate(candidate); 278 } 227 279 228 280 fParticleOutputArray->Add(candidate); … … 235 287 } 236 288 237 vertex = factory->NewCandidate(); 289 nvtx++; 290 238 291 vertex->Position.SetXYZT(vx, vy, dz, dt); 292 293 vertex->ClusterIndex = nvtx; 294 vertex->ClusterNDF = nch; 295 vertex->SumPT2 = sumpt2; 296 vertex->GenSumPT2 = sumpt2; 297 239 298 vertex->IsPU = 1; 240 299 241 300 fVertexOutputArray->Add(vertex); 301 242 302 } 243 303 }
Note:
See TracChangeset
for help on using the changeset viewer.