CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
MixEvtVtxGenerator.cc
Go to the documentation of this file.
1 #ifndef HI_MixEvtVtxGenerator_H
2 #define HI_MixEvtVtxGenerator_H
3 /*
4 * $Date: 2010/11/22 12:41:58 $
5 * $Revision: 1.5 $
6 */
9 
17 
20 
23 
24 #include "TMatrixD.h"
25 
26 using namespace edm;
27 using namespace std;
28 
29 
30 namespace HepMC {
31  class FourVector ;
32 }
33 
34 
36 {
37  public:
38 
39  // ctor & dtor
40  explicit MixEvtVtxGenerator( const edm::ParameterSet& );
41  virtual ~MixEvtVtxGenerator();
42 
43  virtual void produce( edm::Event&, const edm::EventSetup& );
44 
45  virtual HepMC::FourVector* getVertex(edm::Event&);
46  virtual HepMC::FourVector* getRecVertex(edm::Event&);
47 
48  protected:
49 
50  HepMC::FourVector* fVertex ;
51  TMatrixD *boost_;
52 
53  private :
54 
58  std::vector<double> vtxOffset;
59 
60 };
61 
63  : fVertex(0), boost_(0),
64  signalLabel(pset.getParameter<edm::InputTag>("signalLabel")),
65  hiLabel(pset.getParameter<edm::InputTag>("heavyIonLabel")),
66  useRecVertex(pset.exists("useRecVertex")?pset.getParameter<bool>("useRecVertex"):false)
67 
68 {
69  produces<bool>("matchedVertex");
70  vtxOffset.resize(3);
71  if(pset.exists("vtxOffset")) vtxOffset=pset.getParameter< std::vector<double> >("vtxOffset");
72 }
73 
75 {
76  delete fVertex ;
77  if (boost_ != 0 ) delete boost_;
78  // no need since now it's done in HepMCProduct
79  // delete fEvt ;
80 }
81 
82 HepMC::FourVector* MixEvtVtxGenerator::getVertex( Event& evt){
83 
85  evt.getByLabel(hiLabel,input);
86 
87  const HepMC::GenEvent* inev = input->GetEvent();
88  HepMC::GenVertex* genvtx = inev->signal_process_vertex();
89  if(!genvtx){
90  cout<<"No Signal Process Vertex!"<<endl;
91  HepMC::GenEvent::particle_const_iterator pt=inev->particles_begin();
92  HepMC::GenEvent::particle_const_iterator ptend=inev->particles_end();
93  while(!genvtx || ( genvtx->particles_in_size() == 1 && pt != ptend ) ){
94  if(!genvtx) cout<<"No Gen Vertex!"<<endl;
95  if(pt == ptend) cout<<"End reached!"<<endl;
96  genvtx = (*pt)->production_vertex();
97  ++pt;
98  }
99  }
100  double aX,aY,aZ,aT;
101 
102  aX = genvtx->position().x();
103  aY = genvtx->position().y();
104  aZ = genvtx->position().z();
105  aT = genvtx->position().t();
106 
107  if(!fVertex) fVertex = new HepMC::FourVector();
108  fVertex->set(aX,aY,aZ,aT);
109 
110  return fVertex;
111 
112 }
113 
114 
115 HepMC::FourVector* MixEvtVtxGenerator::getRecVertex( Event& evt){
116 
118  evt.getByLabel(hiLabel,input);
119 
120  double aX,aY,aZ;
121 
122  aX = input->begin()->position().x() + vtxOffset[0];
123  aY = input->begin()->position().y() + vtxOffset[1];
124  aZ = input->begin()->position().z() + vtxOffset[2];
125 
126  /*
127  std::cout << "reco::Vertex = " << input->begin()->position().x()
128  << ", " << input->begin()->position().y()
129  << ", " << input->begin()->position().z()
130  << std::endl;
131 
132  std::cout << "offset = " << vtxOffset[0]
133  << ", " << vtxOffset[1]
134  << ", " << vtxOffset[2]
135  << std::endl;
136 
137  std::cout << "embedded GEN vertex = " << aX
138  << ", " << aY << ", " << aZ << std::endl;
139  */
140 
141  if(!fVertex) fVertex = new HepMC::FourVector();
142  fVertex->set(10.0*aX,10.0*aY,10.0*aZ,0.0); // HepMC positions in mm (RECO in cm)
143 
144  return fVertex;
145 
146 }
147 
149 {
150 
151 
152  Handle<HepMCProduct> HepMCEvt ;
153 
154  evt.getByLabel( signalLabel, HepMCEvt ) ;
155 
156  // generate new vertex & apply the shift
157  //
158 
159  HepMCEvt->applyVtxGen( useRecVertex ? getRecVertex(evt) : getVertex(evt) ) ;
160 
161  // HepMCEvt->boostToLab( GetInvLorentzBoost(), "vertex" );
162  // HepMCEvt->boostToLab( GetInvLorentzBoost(), "momentum" );
163 
164  // OK, create a (pseudo)product and put in into edm::Event
165  //
166  auto_ptr<bool> NewProduct(new bool(true)) ;
167  evt.put( NewProduct ,"matchedVertex") ;
168 
169  return ;
170 
171 }
172 
174 
175 #endif
T getParameter(std::string const &) const
virtual HepMC::FourVector * getRecVertex(edm::Event &)
std::vector< double > vtxOffset
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
bool exists(std::string const &parameterName) const
checks if a parameter exists
MixEvtVtxGenerator(const edm::ParameterSet &)
virtual HepMC::FourVector * getVertex(edm::Event &)
edm::InputTag signalLabel
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:85
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:356
HepMC::FourVector * fVertex
virtual void produce(edm::Event &, const edm::EventSetup &)
tuple cout
Definition: gather_cfg.py:121