C++ GenEvent::set_signal_process_vertex方法代码示例

StatusCode ConstPtParticleGun::getNextEvent(HepMC::GenEvent& theEvent) {
  Gaudi::LorentzVector theFourMomentum;
  Gaudi::LorentzVector origin;
  // note: pgdid is set in function generateParticle
  int thePdgId;
  generateParticle(theFourMomentum, origin, thePdgId);

  // create HepMC Vertex --
  // by calling add_vertex(), the hepmc event is given ownership of the vertex
  HepMC::GenVertex* v = new HepMC::GenVertex(HepMC::FourVector(origin.X(), origin.Y(), origin.Z(), origin.T()));
  // create HepMC particle --
  // by calling add_particle_out(), the hepmc vertex is given ownership of the particle
  HepMC::GenParticle* p = new HepMC::GenParticle(
      HepMC::FourVector(theFourMomentum.Px(), theFourMomentum.Py(), theFourMomentum.Pz(), theFourMomentum.E()),
      thePdgId,
      1);  // hepmc status code for final state particle

  v->add_particle_out(p);

  theEvent.add_vertex(v);
  theEvent.set_signal_process_vertex(v);

  return StatusCode::SUCCESS;
}

bool HepMC2_Interface::Sherpa2HepMC(ATOOLS::Blob_List *const blobs,
                                    HepMC::GenEvent& event, double weight)
{
#ifdef USING__HEPMC2__UNITS
  event.use_units(HepMC::Units::GEV,
                  HepMC::Units::MM);
#endif
  event.set_event_number(ATOOLS::rpa->gen.NumberOfGeneratedEvents());
  size_t decid(11);
  std::map<size_t,size_t> decids;
  Blob *sp(blobs->FindFirst(btp::Signal_Process));
  if (sp) {
    Blob_Data_Base *info((*sp)["Decay_Info"]);
    if (info) {
      DecayInfo_Vector decs(info->Get<DecayInfo_Vector>());
      for (size_t i(0);i<decs.size();++i) decids[decs[i]->m_id]=++decid;
    }
  }
  m_blob2genvertex.clear();
  m_particle2genparticle.clear();
  HepMC::GenVertex * vertex;
  std::vector<HepMC::GenParticle*> beamparticles;
  for (ATOOLS::Blob_List::iterator blit=blobs->begin();
       blit!=blobs->end();++blit) {
    if (Sherpa2HepMC(*(blit),vertex,decids)) {
      event.add_vertex(vertex);
      if ((*blit)->Type()==ATOOLS::btp::Signal_Process) {
        if ((**blit)["NLO_subeventlist"]) {
          THROW(fatal_error,"Events containing correlated subtraction events"
                +std::string(" cannot be translated into the full HepMC event")
                +std::string(" format.\n")
                +std::string("   Try 'EVENT_OUTPUT=HepMC_Short' instead."));
        }
        event.set_signal_process_vertex(vertex);
        if((*blit)->NInP()==2) {
          kf_code fl1=(*blit)->InParticle(0)->Flav().HepEvt();
          kf_code fl2=(*blit)->InParticle(1)->Flav().HepEvt();
          double x1=(*blit)->InParticle(0)->Momentum()[0]/rpa->gen.PBeam(0)[0];
          double x2=(*blit)->InParticle(1)->Momentum()[0]/rpa->gen.PBeam(1)[0];
          double q(0.0), p1(0.0), p2(0.0);
          Blob_Data_Base *facscale((**blit)["Factorisation_Scale"]);
	  if (facscale) q=sqrt(facscale->Get<double>());
	  Blob_Data_Base *xf1((**blit)["XF1"]);
          Blob_Data_Base *xf2((**blit)["XF2"]);
          if (xf1) p1=xf1->Get<double>();
          if (xf2) p2=xf2->Get<double>();
	  HepMC::PdfInfo pdfinfo(fl1, fl2, x1, x2, q, p1, p2);
          event.set_pdf_info(pdfinfo);
        }
      }
      else if ((*blit)->Type()==ATOOLS::btp::Beam || 
	       (*blit)->Type()==ATOOLS::btp::Bunch) {
        for (HepMC::GenVertex::particles_in_const_iterator 
	       pit=vertex->particles_in_const_begin();
             pit!=vertex->particles_in_const_end(); ++pit) {
          if ((*pit)->production_vertex()==NULL) {
            beamparticles.push_back(*pit);
          }
        }
      }
    }
  }
  if (beamparticles.size()==2) {
    event.set_beam_particles(beamparticles[0],beamparticles[1]);
  }
  std::vector<double> weights;
  weights.push_back(weight);
  if (sp) {
    Blob_Data_Base *info((*sp)["MEWeight"]);
    if (!info) THROW(fatal_error,"Missing weight info.");
    double meweight(info->Get<double>());
    weights.push_back(meweight);
    Blob_Data_Base *ofni((*sp)["Weight_Norm"]);
    if (!ofni) THROW(fatal_error,"Missing weight normalisation.");
    double weightnorm(ofni->Get<double>());
    weights.push_back(weightnorm);
    ofni=(*sp)["Trials"];
    if (!ofni) THROW(fatal_error,"Missing nof trials.");
    double trials(ofni->Get<double>());
    weights.push_back(trials);
  }
  event.weights()=weights;
  return true;
}