- Timestamp:
- May 26, 2016, 2:54:26 PM (8 years ago)
- Branches:
- ImprovedOutputFile, Timing, dual_readout, llp, master
- Children:
- a889a47
- Parents:
- 971db5b
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
readers/DelphesPythia8.cpp
r971db5b r2acc2e9 152 152 // from pythia8 example 21 153 153 154 void fillParticle(int id, double ee , double thetaIn, double phiIn,154 void fillParticle(int id, double ee_max, double thetaIn, double phiIn, 155 155 Pythia8::Event& event, Pythia8::ParticleData& pdt, Pythia8::Rndm& rndm, bool atRest = false) { 156 156 … … 158 158 event.reset(); 159 159 160 // Select particle mass; where relevant according to Breit-Wigner. 160 // Angles uniform in solid angle. 161 double cThe, sThe, phi, ee; 162 cThe = 2. * rndm.flat() - 1.; 163 sThe = Pythia8::sqrtpos(1. - cThe * cThe); 164 phi = 2. * M_PI * rndm.flat(); 165 ee = pow(10,1+(log10(ee_max)-1)*rndm.flat()); 161 166 double mm = pdt.mSel(id); 162 167 double pp = Pythia8::sqrtpos(ee*ee - mm*mm); 163 168 164 // Special case when particle is supposed to be at rest.165 if (atRest) {166 ee = mm;167 pp = 0.;168 }169 170 // Angles as input or uniform in solid angle.171 double cThe, sThe, phi;172 if (thetaIn >= 0.) {173 cThe = cos(thetaIn);174 sThe = sin(thetaIn);175 phi = phiIn;176 } else {177 cThe = 2. * rndm.flat() - 1.;178 sThe = Pythia8::sqrtpos(1. - cThe * cThe);179 phi = 2. * M_PI * rndm.flat();180 }181 182 169 // Store the particle in the event record. 183 170 event.append( id, 1, 0, 0, pp * sThe * cos(phi), pp * sThe * sin(phi), pp * cThe, ee, mm); … … 185 172 } 186 173 174 void fillPartons(int type, double ee_max, Pythia8::Event& event, Pythia8::ParticleData& pdt, 175 Pythia8::Rndm& rndm) { 176 177 // Reset event record to allow for new event. 178 event.reset(); 179 180 // Angles uniform in solid angle. 181 double cThe, sThe, phi, ee; 182 183 // Information on a q qbar system, to be hadronized. 184 185 cThe = 2. * rndm.flat() - 1.; 186 sThe = Pythia8::sqrtpos(1. - cThe * cThe); 187 phi = 2. * M_PI * rndm.flat(); 188 ee = pow(10,1+(log10(ee_max)-1)*rndm.flat()); 189 double mm = pdt.m0(type); 190 double pp = Pythia8::sqrtpos(ee*ee - mm*mm); 191 if (type == 21) 192 { 193 event.append( 21, 23, 101, 102, pp * sThe * cos(phi), pp * sThe * sin(phi), pp * cThe, ee); 194 event.append( 21, 23, 102, 101, -pp * sThe * cos(phi), -pp * sThe * sin(phi), -pp * cThe, ee); 195 } 196 else 197 { 198 event.append( type, 23, 101, 0, pp * sThe * cos(phi), pp * sThe * sin(phi), pp * cThe, ee, mm); 199 event.append( -type, 23, 0, 101, -pp * sThe * cos(phi), -pp * sThe * sin(phi), -pp * cThe, ee, mm); 200 } 201 } 187 202 188 203 … … 311 326 if (pythia->flag("Main:spareFlag1")) 312 327 { 313 fillParticle( pythia->mode("Main:spareMode1"), 30, -1., 0.,pythia->event, pythia->particleData, pythia->rndm, 0); 328 if (pythia->mode("Main:spareMode1") == 11 || pythia->mode("Main:spareMode1") == 13 || pythia->mode("Main:spareMode1") == 22) 329 { 330 fillParticle( pythia->mode("Main:spareMode1"), pythia->parm("Main:spareParm1"), -1., 0.,pythia->event, pythia->particleData, pythia->rndm, 0); 331 } 332 else fillPartons( pythia->mode("Main:spareMode1"), pythia->parm("Main:spareParm1"), pythia->event, pythia->particleData, pythia->rndm); 314 333 } 315 334
Note:
See TracChangeset
for help on using the changeset viewer.