Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
Merge pull request #17775 from dildick/from-CMSSW_9_1_X_2017-03-04-11…
…00-displaced-muon-gun Add new displaced muon gun
- Loading branch information
Showing
3 changed files
with
240 additions
and
0 deletions.
There are no files selected for viewing
24 changes: 24 additions & 0 deletions
24
Configuration/Generator/python/FlatRandomPtAndDxyGunProducer_cfi.py
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,24 @@ | ||
import FWCore.ParameterSet.Config as cms | ||
|
||
generator = cms.EDProducer("FlatRandomPtAndDxyGunProducer", | ||
PGunParameters = cms.PSet( | ||
PartID = cms.vint32(-13), | ||
MinPt = cms.double(2.00), | ||
MaxPt = cms.double(50.00), | ||
MinEta = cms.double(-3.0), | ||
MaxEta = cms.double(3.0), | ||
MinPhi = cms.double(-3.14159265359), | ||
MaxPhi = cms.double(3.14159265359), | ||
LxyMax = cms.double(3000.0),#make sure most muons generated before Muon system, Gauss distribution | ||
LzMax = cms.double(5000.0),#make sure most muons generated before Muon system, Gauss distribution | ||
ConeRadius = cms.double(1000.0), | ||
ConeH = cms.double(3000.0), | ||
DistanceToAPEX = cms.double(2000.0), | ||
dxyMin = cms.double(0.0), | ||
dxyMax = cms.double(500.0) | ||
), | ||
Verbosity = cms.untracked.int32(0), ## set to 1 (or greater) for printouts | ||
psethack = cms.string('displaced muon'), | ||
AddAntiParticle = cms.bool(True), | ||
firstRun = cms.untracked.uint32(1) | ||
) |
50 changes: 50 additions & 0 deletions
50
IOMC/ParticleGuns/interface/FlatRandomPtAndDxyGunProducer.h
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,50 @@ | ||
#ifndef FlatRandomPtAndDxyGunProducer_H | ||
#define FlatRandomPtAndDxyGunProducer_H | ||
|
||
/** \class FlatRandomPtAndDxyGunProducer | ||
* | ||
* Generates single particle gun in HepMC format | ||
* The d0 is taken by convention equal to dxy | ||
* v1: first version of displaced Muon Gun | ||
* flat in dxy(0, 50), flat in lz(-50, 50) | ||
* v2: bugfixed, unit in HepMC is mm, rather than cm | ||
* v3: check pz*vz>0, make sure muon is not flying | ||
* back to beamspot, Lz is generated by Gauss distribution | ||
* v4: dxy is produced with sign, (-500, 500) [mm] | ||
* v5: check vx*px+vy*py>0, and muons inside beamline-symmetric | ||
* cone is vetoed by lxy^2/c^2<=(vz-d) ^2 && vz-d>0 | ||
* Contact Sven Dildick, Tao Huang | ||
***************************************/ | ||
|
||
#include "IOMC/ParticleGuns/interface/BaseFlatGunProducer.h" | ||
|
||
namespace edm | ||
{ | ||
|
||
class FlatRandomPtAndDxyGunProducer : public BaseFlatGunProducer | ||
{ | ||
|
||
public: | ||
FlatRandomPtAndDxyGunProducer(const ParameterSet & pset); | ||
virtual ~FlatRandomPtAndDxyGunProducer(); | ||
|
||
virtual void produce(Event & e, const EventSetup& es) override; | ||
|
||
private: | ||
|
||
// data members | ||
|
||
double fMinPt ; | ||
double fMaxPt ; | ||
double dxyMin_ ; | ||
double dxyMax_ ; | ||
double lxyMax_ ; | ||
double lzMax_ ; | ||
double ConeRadius_ ; | ||
double ConeH_ ; | ||
double DistanceToAPEX_; | ||
}; | ||
} | ||
|
||
#endif |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,166 @@ | ||
#include <ostream> | ||
|
||
#include "IOMC/ParticleGuns/interface/FlatRandomPtAndDxyGunProducer.h" | ||
|
||
#include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h" | ||
#include "SimDataFormats/GeneratorProducts/interface/GenEventInfoProduct.h" | ||
|
||
#include "FWCore/Framework/interface/Event.h" | ||
#include "FWCore/ParameterSet/interface/ParameterSet.h" | ||
#include "FWCore/ServiceRegistry/interface/Service.h" | ||
#include "FWCore/Utilities/interface/RandomNumberGenerator.h" | ||
|
||
#include "CLHEP/Random/RandFlat.h" | ||
|
||
using namespace edm; | ||
using namespace std; | ||
|
||
FlatRandomPtAndDxyGunProducer::FlatRandomPtAndDxyGunProducer(const ParameterSet& pset) : | ||
BaseFlatGunProducer(pset) | ||
{ | ||
ParameterSet defpset ; | ||
ParameterSet pgun_params = | ||
pset.getParameter<ParameterSet>("PGunParameters") ; | ||
|
||
fMinPt = pgun_params.getParameter<double>("MinPt"); | ||
fMaxPt = pgun_params.getParameter<double>("MaxPt"); | ||
dxyMin_ = pgun_params.getParameter<double>("dxyMin"); | ||
dxyMax_ = pgun_params.getParameter<double>("dxyMax"); | ||
lxyMax_ = pgun_params.getParameter<double>("LxyMax"); | ||
lzMax_ = pgun_params.getParameter<double>("LzMax"); | ||
ConeRadius_ = pgun_params.getParameter<double>("ConeRadius"); | ||
ConeH_ = pgun_params.getParameter<double>("ConeH"); | ||
DistanceToAPEX_ = pgun_params.getParameter<double>("DistanceToAPEX"); | ||
|
||
produces<HepMCProduct>("unsmeared"); | ||
produces<GenEventInfoProduct>(); | ||
} | ||
|
||
FlatRandomPtAndDxyGunProducer::~FlatRandomPtAndDxyGunProducer() | ||
{ | ||
// no need to cleanup GenEvent memory - done in HepMCProduct | ||
} | ||
|
||
void FlatRandomPtAndDxyGunProducer::produce(Event &e, const EventSetup& es) | ||
{ | ||
edm::Service<edm::RandomNumberGenerator> rng; | ||
CLHEP::HepRandomEngine* engine = &rng->getEngine(e.streamID()); | ||
|
||
if ( fVerbosity > 0 ){ | ||
cout << " FlatRandomPtAndDxyGunProducer : Begin New Event Generation" << endl ; | ||
} | ||
// event loop (well, another step in it...) | ||
|
||
// no need to clean up GenEvent memory - done in HepMCProduct | ||
// | ||
|
||
// here re-create fEvt (memory) | ||
// | ||
fEvt = new HepMC::GenEvent() ; | ||
|
||
// now actualy, cook up the event from PDGTable and gun parameters | ||
int barcode = 1 ; | ||
for (unsigned int ip=0; ip<fPartIDs.size(); ++ip){ | ||
|
||
double phi_vtx = 0; | ||
double dxy = 0; | ||
double pt = 0; | ||
double eta = 0; | ||
double px = 0; | ||
double py = 0; | ||
double pz = 0; | ||
double vx = 0; | ||
double vy = 0; | ||
double vz = 0; | ||
double lxy = 0; | ||
|
||
bool passLoop = false; | ||
while (not passLoop) { | ||
|
||
bool passLxy = false; | ||
bool passLz = false; | ||
phi_vtx = CLHEP::RandFlat::shoot(engine, fMinPhi, fMaxPhi); | ||
dxy = CLHEP::RandFlat::shoot(engine, dxyMin_, dxyMax_); | ||
float dxysign = CLHEP::RandFlat::shoot(engine, -1, 1); | ||
if (dxysign < 0) | ||
dxy = -dxy; | ||
pt = CLHEP::RandFlat::shoot(engine, fMinPt, fMaxPt); | ||
px = pt*cos(phi_vtx); | ||
py = pt*sin(phi_vtx); | ||
for (int i=0; i<10000; i++){ | ||
vx = CLHEP::RandFlat::shoot(engine, -lxyMax_, lxyMax_); | ||
vy = (pt*dxy + vx * py)/px; | ||
lxy = sqrt(vx*vx + vy*vy); | ||
if (lxy<abs(lxyMax_) and (vx*px + vy*py)>0){ | ||
passLxy = true; | ||
break; | ||
} | ||
} | ||
eta = CLHEP::RandFlat::shoot(engine, fMinEta, fMaxEta); | ||
pz = pt * sinh(eta); | ||
//vz = fabs(fRandomGaussGenerator->fire(0.0, LzWidth_/2.0)); | ||
float ConeTheta = ConeRadius_/ConeH_; | ||
for (int j=0; j<100; j++){ | ||
vz = CLHEP::RandFlat::shoot(engine, 0.0, lzMax_);// this is abs(vz) | ||
float v0 = vz - DistanceToAPEX_; | ||
if ( v0<=0 or lxy*lxy/(ConeTheta*ConeTheta) > v0*v0 ) { | ||
passLz=true; | ||
break; | ||
} | ||
} | ||
if (pz<0) | ||
vz = -vz; | ||
passLoop = (passLxy and passLz); | ||
|
||
if (passLoop) | ||
break; | ||
} | ||
|
||
HepMC::GenVertex* Vtx1 = new HepMC::GenVertex(HepMC::FourVector(vx,vy,vz)); | ||
|
||
int PartID = fPartIDs[ip] ; | ||
const HepPDT::ParticleData* PData = fPDGTable->particle(HepPDT::ParticleID(abs(PartID))) ; | ||
double mass = PData->mass().value() ; | ||
double energy2= px*px + py*py + pz*pz + mass*mass ; | ||
double energy = sqrt(energy2) ; | ||
HepMC::FourVector p(px,py,pz,energy) ; | ||
HepMC::GenParticle* Part = new HepMC::GenParticle(p,PartID,1); | ||
Part->suggest_barcode( barcode ) ; | ||
barcode++; | ||
Vtx1->add_particle_out(Part); | ||
fEvt->add_vertex(Vtx1) ; | ||
|
||
if ( fAddAntiParticle ) { | ||
HepMC::GenVertex* Vtx2 = new HepMC::GenVertex(HepMC::FourVector(-vx,-vy,-vz)); | ||
HepMC::FourVector ap(-px,-py,-pz,energy) ; | ||
int APartID = -PartID ; | ||
if ( PartID == 22 || PartID == 23 ) | ||
{ | ||
APartID = PartID ; | ||
} | ||
HepMC::GenParticle* APart = new HepMC::GenParticle(ap,APartID,1); | ||
APart->suggest_barcode( barcode ) ; | ||
barcode++ ; | ||
Vtx2->add_particle_out(APart) ; | ||
fEvt->add_vertex(Vtx2) ; | ||
} | ||
} | ||
fEvt->set_event_number(e.id().event()) ; | ||
fEvt->set_signal_process_id(20) ; | ||
|
||
if ( fVerbosity > 0 ){ | ||
fEvt->print() ; | ||
} | ||
|
||
unique_ptr<HepMCProduct> BProduct(new HepMCProduct()) ; | ||
BProduct->addHepMCData( fEvt ); | ||
e.put(std::move(BProduct), "unsmeared"); | ||
|
||
unique_ptr<GenEventInfoProduct> genEventInfo(new GenEventInfoProduct(fEvt)); | ||
e.put(std::move(genEventInfo)); | ||
|
||
if ( fVerbosity > 0 ){ | ||
cout << " FlatRandomPtAndDxyGunProducer : End New Event Generation" << endl ; | ||
fEvt->print(); | ||
} | ||
} |