Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Allow a boost on particles in PythiaMomDauFilter, PythiaFilter and PythiaFilterMultiMother #20140

Merged
merged 7 commits into from Aug 28, 2017
Merged
Show file tree
Hide file tree
Changes from 6 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
16 changes: 16 additions & 0 deletions GeneratorInterface/GenFilters/interface/MCFilterZboostHelper.h
@@ -0,0 +1,16 @@
#ifndef MCFilterZboostHelper_h
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

hi @zc11 - actually I realize CMS has a few implementations of this already. Please use one or adopt the others to use this one - which may require some small refactoring . (https://cmssdt.cern.ch/lxr/search?_filestring=&_string=HepMC%3A%3AFourVector+zboost)

Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for the review.
I've modified the mentioned 3 filters to use MCfilterZboostHelper
Please let me know if it's ok

#define MCFilterZboostHelper_h

#include "HepMC/SimpleVector.h"

namespace HepMC {
class FourVector;
}

namespace MCFilterZboostHelper{

HepMC::FourVector zboost(const HepMC::FourVector&, double);

}

#endif
7 changes: 4 additions & 3 deletions GeneratorInterface/GenFilters/interface/PythiaFilter.h
Expand Up @@ -31,7 +31,6 @@

#include "FWCore/ParameterSet/interface/ParameterSet.h"


//
// class decleration
//
Expand All @@ -42,10 +41,10 @@ namespace edm {
class PythiaFilter : public edm::EDFilter {
public:
explicit PythiaFilter(const edm::ParameterSet&);
~PythiaFilter();
~PythiaFilter() override;


virtual bool filter(edm::Event&, const edm::EventSetup&);
bool filter(edm::Event&, const edm::EventSetup&) override;
private:
// ----------member data ---------------------------

Expand All @@ -67,5 +66,7 @@ class PythiaFilter : public edm::EDFilter {
int status;
int motherID;
int processID;

double betaBoost;
};
#endif
Expand Up @@ -14,7 +14,6 @@

#include "FWCore/ParameterSet/interface/ParameterSet.h"


//
// class decleration
//
Expand All @@ -25,12 +24,12 @@ namespace edm {
class PythiaFilterMultiMother : public edm::EDFilter {
public:
explicit PythiaFilterMultiMother(const edm::ParameterSet&);
~PythiaFilterMultiMother();
~PythiaFilterMultiMother() override;


virtual bool filter(edm::Event&, const edm::EventSetup&);
bool filter(edm::Event&, const edm::EventSetup&) override;
private:
// ----------member data ---------------------------
// ----------member data ---------------------------

edm::EDGetTokenT<edm::HepMCProduct> token_;
int particleID;
Expand All @@ -50,5 +49,7 @@ class PythiaFilterMultiMother : public edm::EDFilter {
int status;
std::vector<int> motherIDs;
int processID;

double betaBoost;
};
#endif
16 changes: 8 additions & 8 deletions GeneratorInterface/GenFilters/interface/PythiaMomDauFilter.h
Expand Up @@ -20,9 +20,7 @@
//
//

namespace edm {
class HepMCProduct;
}

// system include files
#include <memory>

Expand All @@ -35,6 +33,9 @@ namespace edm {

#include "FWCore/ParameterSet/interface/ParameterSet.h"

namespace edm {
class HepMCProduct;
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

same for both of these additions

}

//
// class decleration
Expand All @@ -43,16 +44,14 @@ namespace edm {
class PythiaMomDauFilter : public edm::EDFilter {
public:
explicit PythiaMomDauFilter(const edm::ParameterSet&);
~PythiaMomDauFilter();
~PythiaMomDauFilter() override;


virtual bool filter(edm::Event&, const edm::EventSetup&);
bool filter(edm::Event&, const edm::EventSetup&) override;
private:
// ----------memeber function----------------------

// ----------member data ---------------------------

edm::EDGetTokenT<edm::HepMCProduct> label_;
edm::EDGetTokenT<edm::HepMCProduct> label_;
std::vector<int> dauIDs;
std::vector<int> desIDs;
int particleID;
Expand All @@ -68,6 +67,7 @@ class PythiaMomDauFilter : public edm::EDFilter {
double mom_maxptcut;
double mom_minetacut;
double mom_maxetacut;
double betaBoost;
};
#define PYCOMP pycomp_
extern "C" {
Expand Down
11 changes: 11 additions & 0 deletions GeneratorInterface/GenFilters/src/MCFilterZboostHelper.cc
@@ -0,0 +1,11 @@
#include "GeneratorInterface/GenFilters/interface/MCFilterZboostHelper.h"

HepMC::FourVector MCFilterZboostHelper::zboost(const HepMC::FourVector& mom, double betaBoost) {
//Boost this Lorentz vector (from TLorentzVector::Boost)
double b2 = betaBoost*betaBoost;
double gamma = 1.0 / sqrt(1.0 - b2);
double bp = betaBoost*mom.pz();
double gamma2 = b2 > 0 ? (gamma - 1.0)/b2 : 0.0;

return HepMC::FourVector(mom.px(), mom.py(), mom.pz() + gamma2*bp*betaBoost + gamma*betaBoost*mom.e(), gamma*(mom.e()+bp));
}
17 changes: 9 additions & 8 deletions GeneratorInterface/GenFilters/src/PythiaFilter.cc
@@ -1,5 +1,6 @@

#include "GeneratorInterface/GenFilters/interface/PythiaFilter.h"
#include "GeneratorInterface/GenFilters/interface/MCFilterZboostHelper.h"

#include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
#include <iostream>
Expand All @@ -23,7 +24,8 @@ minphicut(iConfig.getUntrackedParameter("MinPhi", -3.5)),
maxphicut(iConfig.getUntrackedParameter("MaxPhi", 3.5)),
status(iConfig.getUntrackedParameter("Status", 0)),
motherID(iConfig.getUntrackedParameter("MotherID", 0)),
processID(iConfig.getUntrackedParameter("ProcessID", 0))
processID(iConfig.getUntrackedParameter("ProcessID", 0)),
betaBoost(iConfig.getUntrackedParameter("BetaBoost",0.))
{
//now do what ever initialization is needed

Expand Down Expand Up @@ -57,16 +59,16 @@ bool PythiaFilter::filter(edm::Event& iEvent, const edm::EventSetup& iSetup)

for ( HepMC::GenEvent::particle_const_iterator p = myGenEvent->particles_begin();
p != myGenEvent->particles_end(); ++p ) {

rapidity = 0.5*log( ((*p)->momentum().e()+(*p)->momentum().pz()) / ((*p)->momentum().e()-(*p)->momentum().pz()) );
HepMC::FourVector mom = MCFilterZboostHelper::zboost((*p)->momentum(),betaBoost);
rapidity = 0.5*log( (mom.e()+mom.pz()) / (mom.e()-mom.pz()) );

if ( abs((*p)->pdg_id()) == particleID
&& (*p)->momentum().rho() > minpcut
&& (*p)->momentum().rho() < maxpcut
&& mom.rho() > minpcut
&& mom.rho() < maxpcut
&& (*p)->momentum().perp() > minptcut
&& (*p)->momentum().perp() < maxptcut
&& (*p)->momentum().eta() > minetacut
&& (*p)->momentum().eta() < maxetacut
&& mom.eta() > minetacut
&& mom.eta() < maxetacut
&& rapidity > minrapcut
&& rapidity < maxrapcut
&& (*p)->momentum().phi() > minphicut
Expand Down Expand Up @@ -121,4 +123,3 @@ bool PythiaFilter::filter(edm::Event& iEvent, const edm::EventSetup& iSetup)
return true; } else {return false;}

}

16 changes: 9 additions & 7 deletions GeneratorInterface/GenFilters/src/PythiaFilterMultiMother.cc
@@ -1,5 +1,6 @@

#include "GeneratorInterface/GenFilters/interface/PythiaFilterMultiMother.h"
#include "GeneratorInterface/GenFilters/interface/MCFilterZboostHelper.h"

#include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
#include <iostream>
Expand All @@ -23,7 +24,8 @@ minphicut(iConfig.getUntrackedParameter("MinPhi", -3.5)),
maxphicut(iConfig.getUntrackedParameter("MaxPhi", 3.5)),
status(iConfig.getUntrackedParameter("Status", 0)),
motherIDs(iConfig.getUntrackedParameter("MotherIDs", std::vector<int>{0})),
processID(iConfig.getUntrackedParameter("ProcessID", 0))
processID(iConfig.getUntrackedParameter("ProcessID", 0)),
betaBoost(iConfig.getUntrackedParameter("BetaBoost",0.))
{
//now do what ever initialization is needed

Expand Down Expand Up @@ -57,16 +59,16 @@ bool PythiaFilterMultiMother::filter(edm::Event& iEvent, const edm::EventSetup&

for ( HepMC::GenEvent::particle_const_iterator p = myGenEvent->particles_begin();
p != myGenEvent->particles_end(); ++p ) {

rapidity = 0.5*log( ((*p)->momentum().e()+(*p)->momentum().pz()) / ((*p)->momentum().e()-(*p)->momentum().pz()) );
HepMC::FourVector mom = MCFilterZboostHelper::zboost((*p)->momentum(),betaBoost);
rapidity = 0.5*log( (mom.e()+mom.pz()) / (mom.e()-mom.pz()) );

if ( abs((*p)->pdg_id()) == particleID
&& (*p)->momentum().rho() > minpcut
&& (*p)->momentum().rho() < maxpcut
&& mom.rho() > minpcut
&& mom.rho() < maxpcut
&& (*p)->momentum().perp() > minptcut
&& (*p)->momentum().perp() < maxptcut
&& (*p)->momentum().eta() > minetacut
&& (*p)->momentum().eta() < maxetacut
&& mom.eta() > minetacut
&& mom.eta() < maxetacut
&& rapidity > minrapcut
&& rapidity < maxrapcut
&& (*p)->momentum().phi() > minphicut
Expand Down
41 changes: 23 additions & 18 deletions GeneratorInterface/GenFilters/src/PythiaMomDauFilter.cc
@@ -1,6 +1,6 @@

#include "GeneratorInterface/GenFilters/interface/PythiaMomDauFilter.h"

#include "GeneratorInterface/GenFilters/interface/MCFilterZboostHelper.h"

#include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
#include "HepMC/PythiaWrapper6_4.h"
Expand All @@ -12,19 +12,20 @@ using namespace std;

PythiaMomDauFilter::PythiaMomDauFilter(const edm::ParameterSet& iConfig) :
label_(consumes<edm::HepMCProduct>(edm::InputTag(iConfig.getUntrackedParameter("moduleLabel",std::string("generator")),"unsmeared"))),
particleID(iConfig.getUntrackedParameter("ParticleID", 0)),
daughterID(iConfig.getUntrackedParameter("DaughterID", 0)),
chargeconju(iConfig.getUntrackedParameter("ChargeConjugation", true)),
ndaughters(iConfig.getUntrackedParameter("NumberDaughters", 0)),
ndescendants(iConfig.getUntrackedParameter("NumberDescendants", 0)),
minptcut(iConfig.getUntrackedParameter("MinPt", 0.)),
maxptcut(iConfig.getUntrackedParameter("MaxPt", 14000.)),
minetacut(iConfig.getUntrackedParameter("MinEta", -10.)),
maxetacut(iConfig.getUntrackedParameter("MaxEta", 10.)),
mom_minptcut(iConfig.getUntrackedParameter("MomMinPt", 0.)),
mom_maxptcut(iConfig.getUntrackedParameter("MomMaxPt", 14000.)),
mom_minetacut(iConfig.getUntrackedParameter("MomMinEta", -10.)),
mom_maxetacut(iConfig.getUntrackedParameter("MomMaxEta", 10.))
particleID(iConfig.getUntrackedParameter<int>("ParticleID", 0)),
daughterID(iConfig.getUntrackedParameter<int>("DaughterID", 0)),
chargeconju(iConfig.getUntrackedParameter<bool>("ChargeConjugation", true)),
ndaughters(iConfig.getUntrackedParameter<int>("NumberDaughters", 0)),
ndescendants(iConfig.getUntrackedParameter<int>("NumberDescendants", 0)),
minptcut(iConfig.getUntrackedParameter<double>("MinPt", 0.)),
maxptcut(iConfig.getUntrackedParameter<double>("MaxPt", 14000.)),
minetacut(iConfig.getUntrackedParameter<double>("MinEta", -10.)),
maxetacut(iConfig.getUntrackedParameter<double>("MaxEta", 10.)),
mom_minptcut(iConfig.getUntrackedParameter<double>("MomMinPt", 0.)),
mom_maxptcut(iConfig.getUntrackedParameter<double>("MomMaxPt", 14000.)),
mom_minetacut(iConfig.getUntrackedParameter<double>("MomMinEta", -10.)),
mom_maxetacut(iConfig.getUntrackedParameter<double>("MomMaxEta", 10.)),
betaBoost(iConfig.getUntrackedParameter("BetaBoost",0.))
{
//now do what ever initialization is needed
vector<int> defdauID;
Expand Down Expand Up @@ -67,7 +68,8 @@ bool PythiaMomDauFilter::filter(edm::Event& iEvent, const edm::EventSetup& iSetu
for ( HepMC::GenEvent::particle_const_iterator p = myGenEvent->particles_begin(); p != myGenEvent->particles_end(); ++p ) {

if( (*p)->pdg_id() != particleID ) continue ;
if( (*p)->momentum().perp() > mom_minptcut && (*p)->momentum().perp() < mom_maxptcut && (*p)->momentum().eta() > mom_minetacut && (*p)->momentum().eta() < mom_maxetacut ){
HepMC::FourVector mom = MCFilterZboostHelper::zboost((*p)->momentum(),betaBoost);
if( (*p)->momentum().perp() > mom_minptcut && (*p)->momentum().perp() < mom_maxptcut && mom.eta() > mom_minetacut && mom.eta() < mom_maxetacut ){
mom_accepted = true;
ndauac = 0;
ndau = 0;
Expand All @@ -86,7 +88,8 @@ bool PythiaMomDauFilter::filter(edm::Event& iEvent, const edm::EventSetup& iSetu
++ndes;
for( unsigned int i=0; i<desIDs.size(); ++i) {
if( (*des)->pdg_id() != desIDs[i] ) continue ;
if( (*des)->momentum().perp() > minptcut && (*des)->momentum().perp() < maxptcut && (*des)->momentum().eta() > minetacut && (*des)->momentum().eta() < maxetacut ) {
HepMC::FourVector dau_i = MCFilterZboostHelper::zboost((*des)->momentum(),betaBoost);
if( (*des)->momentum().perp() > minptcut && (*des)->momentum().perp() < maxptcut && dau_i.eta() > minetacut && dau_i.eta() < maxetacut ) {
++ndesac;
break;
}
Expand All @@ -110,7 +113,8 @@ bool PythiaMomDauFilter::filter(edm::Event& iEvent, const edm::EventSetup& iSetu
for ( HepMC::GenEvent::particle_const_iterator p = myGenEvent->particles_begin(); p != myGenEvent->particles_end(); ++p ) {

if( (*p)->pdg_id() != -particleID ) continue ;
if( (*p)->momentum().perp() > mom_minptcut && (*p)->momentum().perp() < mom_maxptcut && (*p)->momentum().eta() > mom_minetacut && (*p)->momentum().eta() < mom_maxetacut ){
HepMC::FourVector mom = MCFilterZboostHelper::zboost((*p)->momentum(),betaBoost);
if( (*p)->momentum().perp() > mom_minptcut && (*p)->momentum().perp() < mom_maxptcut && mom.eta() > mom_minetacut && mom.eta() < mom_maxetacut ){
mom_accepted = true;
ndauac = 0;
ndau = 0;
Expand Down Expand Up @@ -141,7 +145,8 @@ bool PythiaMomDauFilter::filter(edm::Event& iEvent, const edm::EventSetup& iSetu
int has_antipart = pydat2.kchg[3-1][pythiaCode-1];
if( has_antipart == 0 ) IDanti = desIDs[i];
if( (*des)->pdg_id() != IDanti ) continue ;
if( (*des)->momentum().perp() > minptcut && (*des)->momentum().perp() < maxptcut && (*des)->momentum().eta() > minetacut && (*des)->momentum().eta() < maxetacut ) {
HepMC::FourVector dau_i = MCFilterZboostHelper::zboost((*des)->momentum(),betaBoost);
if( (*des)->momentum().perp() > minptcut && (*des)->momentum().perp() < maxptcut && dau_i.eta() > minetacut && dau_i.eta() < maxetacut ) {
++ndesac;
break;
}
Expand Down