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

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;
}

//.........这里部分代码省略.........

    for (unsigned int tr = 0; tr < _registeredTriggers.size(); tr++) { 
      bool trigResult = _registeredTriggers[tr]->Apply(event);

      if (verbosity > 2) {
	cout << "PHSartre::process_event trigger: "
	     << _registeredTriggers[tr]->GetName() << "  " << trigResult << endl;
      }

      if (_triggersOR && trigResult) {
	passedTrigger = true;
	break;
      } else if (_triggersAND) {
	andScoreKeeper &= trigResult;
      }
      
      if (verbosity > 2 && !passedTrigger) {
	cout << "PHSartre::process_event - failed trigger: "
	     << _registeredTriggers[tr]->GetName() <<  endl;
      }
    }

    if ((andScoreKeeper && _triggersAND) || (_registeredTriggers.size() == 0)) {
      passedTrigger = true;
    }

  }

  // fill HepMC object with event
  
  HepMC::GenEvent *genevent = new HepMC::GenEvent(HepMC::Units::GEV, HepMC::Units::MM);

  // add some information to the event
  genevent->set_event_number(_eventcount);

  // Set the PDF information
  HepMC::PdfInfo pdfinfo;
  pdfinfo.set_scalePDF(event->Q2);
  genevent->set_pdf_info(pdfinfo); 

  // We would also like to save:
  //
  // event->t;
  // event->x;
  // event->y;
  // event->s;
  // event->W;
  // event->xpom;
  // (event->polarization == transverse ? 0 : 1);
  // (event->diffractiveMode == coherent ? 0 : 1);
  // 
  // but there doesn't seem to be a good place to do so 
  // within the HepMC event information?
  //
  // t, W and Q^2 form a minial set of good variables for diffractive events
  // Maybe what I do is record the input particles to the event at the HepMC
  // vertices and reconstruct the kinematics from there? 
  
  // Create HepMC vertices and add final state particles to them

  // First, the emitter(electron)-virtual photon vertex:

  HepMC::GenVertex* egammavtx = new HepMC::GenVertex(CLHEP::HepLorentzVector(0.0,0.0,0.0,0.0));
  genevent->add_vertex(egammavtx); 

  egammavtx->add_particle_in( 

bool HepMC2_Interface::Sherpa2ShortHepMC(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());
  HepMC::GenVertex * vertex=new HepMC::GenVertex();
  std::vector<HepMC::GenParticle*> beamparticles;
  std::vector<std::pair<HepMC::FourVector,int> > beamparts,
                                                 remnantparts1, remnantparts2;
  Blob *sp(blobs->FindFirst(btp::Signal_Process));
  NLO_subevtlist* subevtlist(NULL);
  ME_wgtinfo* wgtinfo(0);
  
  if (sp) {
    Blob_Data_Base* seinfo=(*sp)["ME_wgtinfo"];
    if (seinfo) wgtinfo=seinfo->Get<ME_wgtinfo*>();
   
    Blob_Data_Base * bdb((*sp)["NLO_subeventlist"]);
    if (bdb) subevtlist=bdb->Get<NLO_subevtlist*>();
  }
  for (ATOOLS::Blob_List::iterator blit=blobs->begin();
       blit!=blobs->end();++blit) {
    Blob* blob=*blit;
    for (int i=0;i<blob->NInP();i++) {
      if (blob->InParticle(i)->ProductionBlob()==NULL) {
        Particle* parton=blob->InParticle(i);
        ATOOLS::Vec4D mom  = parton->Momentum();
        HepMC::FourVector momentum(mom[1],mom[2],mom[3],mom[0]);
        HepMC::GenParticle* inpart = 
	  new HepMC::GenParticle(momentum,parton->Flav().HepEvt(),2);
        vertex->add_particle_in(inpart);
        // distinct because SHRIMPS has no bunches for some reason
        if (blob->Type()==btp::Beam || blob->Type()==btp::Bunch) {
          beamparticles.push_back(inpart);
          beamparts.push_back(std::make_pair(momentum,parton->Flav().HepEvt()));
        }
      }
    }
    for (int i=0;i<blob->NOutP();i++) {
      if (blob->OutParticle(i)->DecayBlob()==NULL) {
        Particle* parton=blob->OutParticle(i);
        ATOOLS::Vec4D mom  = parton->Momentum();
        HepMC::FourVector momentum(mom[1],mom[2],mom[3],mom[0]);
        HepMC::GenParticle* outpart = 
	  new HepMC::GenParticle(momentum,parton->Flav().HepEvt(),1);
        vertex->add_particle_out(outpart);
        if (blob->Type()==btp::Beam) {
          if      (mom[3]>0) 
	    remnantparts1.push_back(std::make_pair(momentum,
						   parton->Flav().HepEvt()));
          else if (mom[3]<0) 
	    remnantparts2.push_back(std::make_pair(momentum,
						   parton->Flav().HepEvt()));
          else THROW(fatal_error,"Ill defined beam remnants.");
        }
      }
    }
    
    if ((*blit)->Type()==ATOOLS::btp::Signal_Process) {
      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);
      }
    }
  }
  event.add_vertex(vertex);
  if (beamparticles.size()==2) {
    event.set_beam_particles(beamparticles[0],beamparticles[1]);
  }
  std::vector<double> weights; weights.push_back(weight);
  if (sp && !subevtlist) {
    Blob_Data_Base *info;
    info=((*sp)["MEWeight"]);
    if (!info) THROW(fatal_error,"Missing weight info.");
    double meweight(info->Get<double>());
    weights.push_back(meweight);
    info=((*sp)["Weight_Norm"]);
    if (!info) THROW(fatal_error,"Missing weight normalisation.");
    double weightnorm(info->Get<double>());
    weights.push_back(weightnorm);
    info=(*sp)["Trials"];
    if (!info) THROW(fatal_error,"Missing nof trials.");
    double trials(info->Get<double>());
    weights.push_back(trials);
    //alphaS value && power
    double rscale2 = (*sp)["Renormalization_Scale"]->Get<double>();
//.........这里部分代码省略.........