[u/mrichter/AliRoot.git] / THerwig / THerwig6.cxx


// definition of c++ Class THerwig6 to be used in ROOT
// this is a c++ interface to the F77 Herwig6 program
// author: j. g. contreras jgcn@moni.mda.cinvestav.mx
// date: december 22, 2000

#include <stream.h>
#include "THerwig6.h"
#include "TClonesArray.h"
#include "TParticle.h"


ClassImp(THerwig6)

extern "C" {
  void*  herwig6_common_block_address_(char*, int len);
  void   herwig6_open_fortran_file_ (int* lun, char* name, int);
  void   herwig6_close_fortran_file_(int* lun);
}


THerwig6::THerwig6() : TGenerator("Herwig6","Herwig6") {

// THerwig6 constructor: creates a TClonesArray in which it will store all
// particles. Note that there may be only one functional THerwig6 object
// at a time, so it's not use to create more than one instance of it.

  delete fParticles; // was allocated as TObjArray in TGenerator

  fParticles = new TClonesArray("TParticle",50);

  // initialize common-blocks
  fHepevt = (Hepevt_t*) herwig6_common_block_address_((char*)"HEPEVT",6);
  fHwbeam = (Hwbeam_t*) herwig6_common_block_address_((char*)"HWBEAM",6);
  fHwbmch = (Hwbmch_t*) herwig6_common_block_address_((char*)"HWBMCH",6);
  fHwproc = (Hwproc_t*) herwig6_common_block_address_((char*)"HWPROC",6);
  fHwpram = (Hwpram_t*) herwig6_common_block_address_((char*)"HWPRAM",6);
  fHwprch = (Hwprch_t*) herwig6_common_block_address_((char*)"HWPRCH",6);
  fHwpart = (Hwpart_t*) herwig6_common_block_address_((char*)"HWPART",6);
  fHwparp = (Hwparp_t*) herwig6_common_block_address_((char*)"HWPARP",6);
  fHwbosc = (Hwbosc_t*) herwig6_common_block_address_((char*)"HWBOSC",6);
  fHwparc = (Hwparc_t*) herwig6_common_block_address_((char*)"HWPARC",6);
  fHwbrch = (Hwbrch_t*) herwig6_common_block_address_((char*)"HWBRCH",6);
  fHwevnt = (Hwevnt_t*) herwig6_common_block_address_((char*)"HWEVNT",6);
  fHwhard = (Hwhard_t*) herwig6_common_block_address_((char*)"HWHARD",6);
  fHwprop = (Hwprop_t*) herwig6_common_block_address_((char*)"HWPROP",6);
  fHwunam = (Hwunam_t*) herwig6_common_block_address_((char*)"HWUNAM",6);
  fHwupdt = (Hwupdt_t*) herwig6_common_block_address_((char*)"HWUPDT",6);
  fHwuwts = (Hwuwts_t*) herwig6_common_block_address_((char*)"HWUWTS",6);
  fHwuclu = (Hwuclu_t*) herwig6_common_block_address_((char*)"HWUCLU",6);
  fHwdist = (Hwdist_t*) herwig6_common_block_address_((char*)"HWDIST",6);
  fHwqdks = (Hwqdks_t*) herwig6_common_block_address_((char*)"HWQDKS",6);
  fHwusud = (Hwusud_t*) herwig6_common_block_address_((char*)"HWUSUD",6);
  fHwsusy = (Hwsusy_t*) herwig6_common_block_address_((char*)"HWSUSY",6);
  fHwrpar = (Hwrpar_t*) herwig6_common_block_address_((char*)"HWRPAR",6);
  fHwminb = (Hwminb_t*) herwig6_common_block_address_((char*)"HWMINB",6);
  fHwclus = (Hwclus_t*) herwig6_common_block_address_((char*)"HWCLUS",6);
 }

//------------------------------------------------------------------------------
 THerwig6::~THerwig6()
 {
// Destroys the object, deletes and disposes all TParticles currently on list.

    if (fParticles) {
      fParticles->Delete();
      delete fParticles;
      fParticles = 0;
   }
 }

//______________________________________________________________________________
void THerwig6::GenerateEvent() 
{
  
  //initialize event
  hwuine_();
  // generate hard subprocess
  hwepro_();
  // generate parton cascades
  hwbgen_();
  // do heavy objects decay
  hwdhob_();
  // do cluster formation
  hwcfor_();
  // do cluster decays
  hwcdec_();
  // do unstable particle decays
  hwdhad_();
  // do heavy flavor hadrons decay
  hwdhvy_();
  // add soft underlying event
  hwmevt_();
  // finish event
  hwufne_();
}

//______________________________________________________________________________
void THerwig6::OpenFortranFile(int lun, char* name) {
  herwig6_open_fortran_file_(&lun, name, strlen(name));
}

//______________________________________________________________________________
void THerwig6::CloseFortranFile(int lun) {
  herwig6_close_fortran_file_(&lun);
}

void THerwig6::Initialize(const char *beam, const char *target, double pbeam1, double pbeam2, int iproc)

{
  // perform the initialization for Herwig6
  // sets correct title. 
  // after calling this method all parameters are set to their default
  // values. If you want to modify any parameter you have to set the new
  // value after calling Initialize and before PrepareRun.

   char  cbeam[8];
   strncpy(cbeam,beam,8);
   char  ctarget[8];
   strncpy(ctarget,target,8);
   printf("\n Initializing Herwig !! \n");
   if ( (!strncmp(beam, "E+"    ,2)) &&
        (!strncmp(beam, "E-"    ,2)) &&
        (!strncmp(beam, "MU+"   ,3)) &&
        (!strncmp(beam, "MU-"   ,3)) &&
        (!strncmp(beam, "NUE"   ,3)) &&
        (!strncmp(beam, "NUEB"  ,4)) &&
        (!strncmp(beam, "NUMU"  ,4)) &&
        (!strncmp(beam, "NMUB"  ,4)) &&
        (!strncmp(beam, "NTAU"  ,4)) &&
        (!strncmp(beam, "NTAB"  ,4)) &&
        (!strncmp(beam, "GAMA"  ,4)) &&
        (!strncmp(beam, "P       ",8)) &&
        (!strncmp(beam, "PBAR    ",8)) &&
        (!strncmp(beam, "N"     ,1)) &&
        (!strncmp(beam, "NBAR"  ,4)) &&
        (!strncmp(beam, "PI+"   ,3)) &&
        (!strncmp(beam, "PI-"   ,3)) ) {
      printf("WARNING! In THerwig6:Initialize():\n");
      printf(" specified beam=%s is unrecognized .\n",beam);
      printf(" resetting to \"P\" .");
      sprintf(cbeam,"P");
   }

   if ( (!strncmp(target, "E+"    ,2)) &&
        (!strncmp(target, "E-"    ,2)) &&
        (!strncmp(target, "MU+"   ,3)) &&
        (!strncmp(target, "MU-"   ,3)) &&
        (!strncmp(target, "NUE"   ,3)) &&
        (!strncmp(target, "NUEB"  ,4)) &&
        (!strncmp(target, "NUMU"  ,4)) &&
        (!strncmp(target, "NMUB"  ,4)) &&
        (!strncmp(target, "NTAU"  ,4)) &&
        (!strncmp(target, "NTAB"  ,4)) &&
        (!strncmp(target, "GAMA"  ,4)) &&
        (!strncmp(target, "P       ",8)) &&
        (!strncmp(target, "PBAR    ",8)) &&
        (!strncmp(target, "N"     ,1)) &&
        (!strncmp(target, "NBAR"  ,4)) &&
        (!strncmp(target, "PI+"   ,3)) &&
        (!strncmp(target, "PI-"   ,3)) ) {
      printf("WARNING! In THerwig6:Initialize():\n");
      printf(" specified target=%s is unrecognized .\n",target);
      printf(" resetting to \"P\" .");
      sprintf(ctarget,"P");
   }

   // initialization:
   // type of beams
   strncpy(fHwbmch->PART1,beam,8);
   strncpy(fHwbmch->PART2,target,8);
   // momentum of beams
   fHwproc->PBEAM1=pbeam1;
   fHwproc->PBEAM2=pbeam2;
   // process to generate
   fHwproc->IPROC=iproc;
   // not used in the class definition 
   fHwproc->MAXEV=1;

   // reset all parameters
   hwigin_();
   
   // set correct title
   char atitle[132];
   double win=pbeam1+pbeam2;
   printf("\n %s - %s at %g GeV",beam,target,win);
   sprintf(atitle,"%s-%s at %g GeV",cbeam,ctarget,win);
   SetTitle(atitle); 
}

void THerwig6::PrepareRun()
{
  // compute parameter dependent constants
  hwuinc_();
  // initialize elementary processes
  hweini_();
}

//______________________________________________________________________________
TObjArray* THerwig6::ImportParticles(Option_t *option)
{
//
//  Default primary creation method. It reads the /HEPEVT/ common block which
//  has been filled by the GenerateEvent method. If the event generator does
//  not use the HEPEVT common block, This routine has to be overloaded by
//  the subclasses.
//  The default action is to store only the stable particles (ISTHEP = 1)
//  This can be demanded explicitly by setting the option = "Final"
//  If the option = "All", all the particles are stored.
//
  fParticles->Clear();
  Int_t numpart = fHepevt->NHEP;
  TClonesArray &a = *((TClonesArray*)fParticles);
  if (!strcmp(option,"") || !strcmp(option,"Final")) {
    for (Int_t i = 0; i<=numpart; i++) {
      if (fHepevt->ISTHEP[i] == 1) {
//
//  Use the common block values for the TParticle constructor
//
        new(a[i]) TParticle(
                                   fHepevt->IDHEP[i],
                                   fHepevt->ISTHEP[i],
                                   fHepevt->JMOHEP[i][0]-1,
                                   fHepevt->JMOHEP[i][1]-1,
                                   fHepevt->JDAHEP[i][0]-1,
                                   fHepevt->JDAHEP[i][1]-1,

                                   fHepevt->PHEP[i][0],
                                   fHepevt->PHEP[i][1],
                                   fHepevt->PHEP[i][2],
                                   fHepevt->PHEP[i][3],
                                   fHepevt->VHEP[i][0],
                                   fHepevt->VHEP[i][1],
                                   fHepevt->VHEP[i][2],
                                   fHepevt->VHEP[i][3]);
        }
     }
  }
  else if (!strcmp(option,"All")) {
    for (Int_t i = 0; i<=numpart; i++) {
      new(a[i]) TParticle(
                                   fHepevt->IDHEP[i],
                                   fHepevt->ISTHEP[i],
                                   fHepevt->JMOHEP[i][0]-1,
                                   fHepevt->JMOHEP[i][1]-1,
                                   fHepevt->JDAHEP[i][0]-1,
                                   fHepevt->JDAHEP[i][1]-1,

                                   fHepevt->PHEP[i][0],
                                   fHepevt->PHEP[i][1],
                                   fHepevt->PHEP[i][2],
                                   fHepevt->PHEP[i][3],
                                   fHepevt->VHEP[i][0],
                                   fHepevt->VHEP[i][1],
                                   fHepevt->VHEP[i][2],
                                   fHepevt->VHEP[i][3]);
    }
  }
  return fParticles;
}

//______________________________________________________________________________
Int_t THerwig6::ImportParticles(TClonesArray *particles, Option_t *option)
{
//
//  Default primary creation method. It reads the /HEPEVT/ common block which
//  has been filled by the GenerateEvent method. If the event generator does
//  not use the HEPEVT common block, This routine has to be overloaded by
//  the subclasses.
//  The function loops on the generated particles and store them in
//  the TClonesArray pointed by the argument particles.
//  The default action is to store only the stable particles (ISTHEP = 1)
//  This can be demanded explicitly by setting the option = "Final"
//  If the option = "All", all the particles are stored.
//
  if (particles == 0) return 0;
  TClonesArray &Particles = *particles;
  Particles.Clear();
  Int_t numpart = fHepevt->NHEP;
  if (!strcmp(option,"") || !strcmp(option,"Final")) {
    for (Int_t i = 0; i< numpart; i++) {
      if (fHepevt->ISTHEP[i] == 1) {
//
//  Use the common block values for the TParticle constructor
//
        new(Particles[i]) TParticle(
                                   fHepevt->IDHEP[i],
                                   fHepevt->ISTHEP[i],
                                   fHepevt->JMOHEP[i][0]-1,
                                   fHepevt->JMOHEP[i][1]-1,
                                   fHepevt->JDAHEP[i][0]-1,
                                   fHepevt->JDAHEP[i][1]-1,

                                   fHepevt->PHEP[i][0],
                                   fHepevt->PHEP[i][1],
                                   fHepevt->PHEP[i][2],
                                   fHepevt->PHEP[i][3],
                                   fHepevt->VHEP[i][0],
                                   fHepevt->VHEP[i][1],
                                   fHepevt->VHEP[i][2],
                                   fHepevt->VHEP[i][3]);
        }
     }
  }
  else if (!strcmp(option,"All")) {
    for (Int_t i = 0; i< numpart; i++) {
      new(Particles[i]) TParticle(
                                   fHepevt->IDHEP[i],
                                   fHepevt->ISTHEP[i],
                                   fHepevt->JMOHEP[i][0]-1,
                                   fHepevt->JMOHEP[i][1]-1,
                                   fHepevt->JDAHEP[i][0]-1,
                                   fHepevt->JDAHEP[i][1]-1,

                                   fHepevt->PHEP[i][0],
                                   fHepevt->PHEP[i][1],
                                   fHepevt->PHEP[i][2],
                                   fHepevt->PHEP[i][3],
                                   fHepevt->VHEP[i][0],
                                   fHepevt->VHEP[i][1],
                                   fHepevt->VHEP[i][2],
                                   fHepevt->VHEP[i][3]);
    }
  }
  return numpart;
}

void THerwig6::Hwigin()
{
  hwigin_();
}

void THerwig6::Hwuinc()
{
  hwuinc_();
}

void THerwig6::Hwusta(char* name)

{
  hwusta_(name,8);
}

void THerwig6::Hweini()

{
  hweini_();
}

void THerwig6::Hwuine()

{
  hwuine_();
}

void THerwig6::Hwepro()

{
  hwepro_();
}

void THerwig6::Hwbgen()

{
  hwbgen_();
}

void THerwig6::Hwdhob()

{
  hwdhob_();
}

void THerwig6::Hwcfor()

{
  hwcfor_();
}

void THerwig6::Hwcdec()

{
  hwcdec_();
}

void THerwig6::Hwdhad()

{
  hwdhad_();
}

void THerwig6::Hwdhvy()

{
  hwdhvy_();
}

void THerwig6::Hwmevt()

{
  hwmevt_();
}

void THerwig6::Hwufne()

{
  hwufne_();
}

void THerwig6::Hwefin()

{
  hwefin_();
}

void THerwig6::SetupTest()
{
  // exampe of running herwig and generating one event
  // after changing some options
  Initialize("P","PBAR",900.,900.,1500);
  // here you can set some parameters
  SetPTMIN(15.); // Min pt in hadronic jet production
  SetYJMIN(-4.); // Min jet rapidity
  SetYJMAX(4.);  // Max jet rapidity
  // after you set your wished parameters
  // herwig can do its work
  PrepareRun();
  int n_ev_to_generate=1;
  for (int i=0;i<n_ev_to_generate;i++)
    {
      GenerateEvent();
      // do your stuff. For ex:
      int n_of_par=GetNumberOfParticles(); // from TGenerator
      for (int j=0; j<n_of_par; j++)
	{
	  TParticle* p=GetParticle(j);
	  // here you do whatever you want with the particle
	  p->Print();
	};
    };
}
Commit	Line	Data
2d8d4543	1
	2	// definition of c++ Class THerwig6 to be used in ROOT
	3	// this is a c++ interface to the F77 Herwig6 program
	4	// author: j. g. contreras jgcn@moni.mda.cinvestav.mx
	5	// date: december 22, 2000
	6
	7	#include <stream.h>
	8	#include "THerwig6.h"
	9	#include "TClonesArray.h"
	10	#include "TParticle.h"
	11
	12
	13	ClassImp(THerwig6)
	14
	15	extern "C" {
	16	void* herwig6_common_block_address_(char*, int len);
	17	void herwig6_open_fortran_file_ (int* lun, char* name, int);
	18	void herwig6_close_fortran_file_(int* lun);
	19	}
	20
	21
	22	THerwig6::THerwig6() : TGenerator("Herwig6","Herwig6") {
	23
	24	// THerwig6 constructor: creates a TClonesArray in which it will store all
	25	// particles. Note that there may be only one functional THerwig6 object
	26	// at a time, so it's not use to create more than one instance of it.
	27
	28	delete fParticles; // was allocated as TObjArray in TGenerator
	29
	30	fParticles = new TClonesArray("TParticle",50);
	31
	32	// initialize common-blocks
	33	fHepevt = (Hepevt_t) herwig6_common_block_address_((char)"HEPEVT",6);
	34	fHwbeam = (Hwbeam_t) herwig6_common_block_address_((char)"HWBEAM",6);
	35	fHwbmch = (Hwbmch_t) herwig6_common_block_address_((char)"HWBMCH",6);
	36	fHwproc = (Hwproc_t) herwig6_common_block_address_((char)"HWPROC",6);
	37	fHwpram = (Hwpram_t) herwig6_common_block_address_((char)"HWPRAM",6);
	38	fHwprch = (Hwprch_t) herwig6_common_block_address_((char)"HWPRCH",6);
	39	fHwpart = (Hwpart_t) herwig6_common_block_address_((char)"HWPART",6);
	40	fHwparp = (Hwparp_t) herwig6_common_block_address_((char)"HWPARP",6);
	41	fHwbosc = (Hwbosc_t) herwig6_common_block_address_((char)"HWBOSC",6);
	42	fHwparc = (Hwparc_t) herwig6_common_block_address_((char)"HWPARC",6);
	43	fHwbrch = (Hwbrch_t) herwig6_common_block_address_((char)"HWBRCH",6);
	44	fHwevnt = (Hwevnt_t) herwig6_common_block_address_((char)"HWEVNT",6);
	45	fHwhard = (Hwhard_t) herwig6_common_block_address_((char)"HWHARD",6);
	46	fHwprop = (Hwprop_t) herwig6_common_block_address_((char)"HWPROP",6);
	47	fHwunam = (Hwunam_t) herwig6_common_block_address_((char)"HWUNAM",6);
	48	fHwupdt = (Hwupdt_t) herwig6_common_block_address_((char)"HWUPDT",6);
	49	fHwuwts = (Hwuwts_t) herwig6_common_block_address_((char)"HWUWTS",6);
	50	fHwuclu = (Hwuclu_t) herwig6_common_block_address_((char)"HWUCLU",6);
	51	fHwdist = (Hwdist_t) herwig6_common_block_address_((char)"HWDIST",6);
	52	fHwqdks = (Hwqdks_t) herwig6_common_block_address_((char)"HWQDKS",6);
	53	fHwusud = (Hwusud_t) herwig6_common_block_address_((char)"HWUSUD",6);
	54	fHwsusy = (Hwsusy_t) herwig6_common_block_address_((char)"HWSUSY",6);
	55	fHwrpar = (Hwrpar_t) herwig6_common_block_address_((char)"HWRPAR",6);
	56	fHwminb = (Hwminb_t) herwig6_common_block_address_((char)"HWMINB",6);
	57	fHwclus = (Hwclus_t) herwig6_common_block_address_((char)"HWCLUS",6);
	58	}
	59
	60	//------------------------------------------------------------------------------
	61	THerwig6::~THerwig6()
	62	{
	63	// Destroys the object, deletes and disposes all TParticles currently on list.
	64
65	if (fParticles) {
66	fParticles->Delete();
67	delete fParticles;
68	fParticles = 0;
69	}
70	}
71
72	//______________________________________________________________________________
73	void THerwig6::GenerateEvent()
74	{
75
76	//initialize event
77	hwuine_();
78	// generate hard subprocess
79	hwepro_();
80	// generate parton cascades
81	hwbgen_();
82	// do heavy objects decay
83	hwdhob_();
84	// do cluster formation
85	hwcfor_();
86	// do cluster decays
87	hwcdec_();
88	// do unstable particle decays
89	hwdhad_();
90	// do heavy flavor hadrons decay
91	hwdhvy_();
92	// add soft underlying event
93	hwmevt_();
94	// finish event
95	hwufne_();
96	}
97
98	//______________________________________________________________________________
99	void THerwig6::OpenFortranFile(int lun, char* name) {
100	herwig6_open_fortran_file_(&lun, name, strlen(name));
101	}
102
103	//______________________________________________________________________________
104	void THerwig6::CloseFortranFile(int lun) {
105	herwig6_close_fortran_file_(&lun);
106	}
107
108	void THerwig6::Initialize(const char beam, const char target, double pbeam1, double pbeam2, int iproc)
109
110	{
111	// perform the initialization for Herwig6
112	// sets correct title.
113	// after calling this method all parameters are set to their default
114	// values. If you want to modify any parameter you have to set the new
115	// value after calling Initialize and before PrepareRun.
116
117	char cbeam[8];
118	strncpy(cbeam,beam,8);
119	char ctarget[8];
120	strncpy(ctarget,target,8);
121	printf("\n Initializing Herwig !! \n");
122	if ( (!strncmp(beam, "E+" ,2)) &&
123	(!strncmp(beam, "E-" ,2)) &&
124	(!strncmp(beam, "MU+" ,3)) &&
125	(!strncmp(beam, "MU-" ,3)) &&
126	(!strncmp(beam, "NUE" ,3)) &&
127	(!strncmp(beam, "NUEB" ,4)) &&
128	(!strncmp(beam, "NUMU" ,4)) &&
129	(!strncmp(beam, "NMUB" ,4)) &&
130	(!strncmp(beam, "NTAU" ,4)) &&
131	(!strncmp(beam, "NTAB" ,4)) &&
132	(!strncmp(beam, "GAMA" ,4)) &&
133	(!strncmp(beam, "P ",8)) &&
134	(!strncmp(beam, "PBAR ",8)) &&
135	(!strncmp(beam, "N" ,1)) &&
136	(!strncmp(beam, "NBAR" ,4)) &&
137	(!strncmp(beam, "PI+" ,3)) &&
138	(!strncmp(beam, "PI-" ,3)) ) {
139	printf("WARNING! In THerwig6:Initialize():\n");
140	printf(" specified beam=%s is unrecognized .\n",beam);
141	printf(" resetting to \"P\" .");
142	sprintf(cbeam,"P");
143	}
144
145	if ( (!strncmp(target, "E+" ,2)) &&
146	(!strncmp(target, "E-" ,2)) &&
147	(!strncmp(target, "MU+" ,3)) &&
148	(!strncmp(target, "MU-" ,3)) &&
149	(!strncmp(target, "NUE" ,3)) &&
150	(!strncmp(target, "NUEB" ,4)) &&
151	(!strncmp(target, "NUMU" ,4)) &&
152	(!strncmp(target, "NMUB" ,4)) &&
153	(!strncmp(target, "NTAU" ,4)) &&
154	(!strncmp(target, "NTAB" ,4)) &&
155	(!strncmp(target, "GAMA" ,4)) &&
156	(!strncmp(target, "P ",8)) &&
157	(!strncmp(target, "PBAR ",8)) &&
158	(!strncmp(target, "N" ,1)) &&
159	(!strncmp(target, "NBAR" ,4)) &&
160	(!strncmp(target, "PI+" ,3)) &&
161	(!strncmp(target, "PI-" ,3)) ) {
162	printf("WARNING! In THerwig6:Initialize():\n");
163	printf(" specified target=%s is unrecognized .\n",target);
164	printf(" resetting to \"P\" .");
165	sprintf(ctarget,"P");
166	}
167
168	// initialization:
169	// type of beams
170	strncpy(fHwbmch->PART1,beam,8);
171	strncpy(fHwbmch->PART2,target,8);
172	// momentum of beams
173	fHwproc->PBEAM1=pbeam1;
174	fHwproc->PBEAM2=pbeam2;
175	// process to generate
176	fHwproc->IPROC=iproc;
177	// not used in the class definition
178	fHwproc->MAXEV=1;
179
180	// reset all parameters
181	hwigin_();
182
183	// set correct title
184	char atitle[132];
185	double win=pbeam1+pbeam2;
186	printf("\n %s - %s at %g GeV",beam,target,win);
187	sprintf(atitle,"%s-%s at %g GeV",cbeam,ctarget,win);
188	SetTitle(atitle);
189	}
190
191	void THerwig6::PrepareRun()
192	{
193	// compute parameter dependent constants
194	hwuinc_();
195	// initialize elementary processes
196	hweini_();
197	}
198
199	//______________________________________________________________________________
200	TObjArray* THerwig6::ImportParticles(Option_t *option)
201	{
202	//
203	// Default primary creation method. It reads the /HEPEVT/ common block which
204	// has been filled by the GenerateEvent method. If the event generator does
205	// not use the HEPEVT common block, This routine has to be overloaded by
206	// the subclasses.
207	// The default action is to store only the stable particles (ISTHEP = 1)
208	// This can be demanded explicitly by setting the option = "Final"
209	// If the option = "All", all the particles are stored.
210	//
211	fParticles->Clear();
212	Int_t numpart = fHepevt->NHEP;
213	TClonesArray &a = ((TClonesArray)fParticles);
214	if (!strcmp(option,"") \|\| !strcmp(option,"Final")) {
215	for (Int_t i = 0; i<=numpart; i++) {
216	if (fHepevt->ISTHEP[i] == 1) {
217	//
218	// Use the common block values for the TParticle constructor
219	//
220	new(a[i]) TParticle(
221	fHepevt->IDHEP[i],
222	fHepevt->ISTHEP[i],
223	fHepevt->JMOHEP[i][0]-1,
224	fHepevt->JMOHEP[i][1]-1,
225	fHepevt->JDAHEP[i][0]-1,
226	fHepevt->JDAHEP[i][1]-1,
227
228	fHepevt->PHEP[i][0],
229	fHepevt->PHEP[i][1],
230	fHepevt->PHEP[i][2],
231	fHepevt->PHEP[i][3],
232	fHepevt->VHEP[i][0],
233	fHepevt->VHEP[i][1],
234	fHepevt->VHEP[i][2],
235	fHepevt->VHEP[i][3]);
236	}
237	}
238	}
239	else if (!strcmp(option,"All")) {
240	for (Int_t i = 0; i<=numpart; i++) {
241	new(a[i]) TParticle(
242	fHepevt->IDHEP[i],
243	fHepevt->ISTHEP[i],
244	fHepevt->JMOHEP[i][0]-1,
245	fHepevt->JMOHEP[i][1]-1,
246	fHepevt->JDAHEP[i][0]-1,
247	fHepevt->JDAHEP[i][1]-1,
248
249	fHepevt->PHEP[i][0],
250	fHepevt->PHEP[i][1],
251	fHepevt->PHEP[i][2],
252	fHepevt->PHEP[i][3],
253	fHepevt->VHEP[i][0],
254	fHepevt->VHEP[i][1],
255	fHepevt->VHEP[i][2],
256	fHepevt->VHEP[i][3]);
257	}
258	}
259	return fParticles;
260	}
261
262	//______________________________________________________________________________
263	Int_t THerwig6::ImportParticles(TClonesArray particles, Option_t option)
264	{
265	//
266	// Default primary creation method. It reads the /HEPEVT/ common block which
267	// has been filled by the GenerateEvent method. If the event generator does
268	// not use the HEPEVT common block, This routine has to be overloaded by
269	// the subclasses.
270	// The function loops on the generated particles and store them in
271	// the TClonesArray pointed by the argument particles.
272	// The default action is to store only the stable particles (ISTHEP = 1)
273	// This can be demanded explicitly by setting the option = "Final"
274	// If the option = "All", all the particles are stored.
275	//
276	if (particles == 0) return 0;
277	TClonesArray &Particles = *particles;
278	Particles.Clear();
279	Int_t numpart = fHepevt->NHEP;
280	if (!strcmp(option,"") \|\| !strcmp(option,"Final")) {
281	for (Int_t i = 0; i< numpart; i++) {
282	if (fHepevt->ISTHEP[i] == 1) {
283	//
284	// Use the common block values for the TParticle constructor
285	//
286	new(Particles[i]) TParticle(
287	fHepevt->IDHEP[i],
288	fHepevt->ISTHEP[i],
289	fHepevt->JMOHEP[i][0]-1,
290	fHepevt->JMOHEP[i][1]-1,
291	fHepevt->JDAHEP[i][0]-1,
292	fHepevt->JDAHEP[i][1]-1,
293
294	fHepevt->PHEP[i][0],
295	fHepevt->PHEP[i][1],
296	fHepevt->PHEP[i][2],
297	fHepevt->PHEP[i][3],
298	fHepevt->VHEP[i][0],
299	fHepevt->VHEP[i][1],
300	fHepevt->VHEP[i][2],
301	fHepevt->VHEP[i][3]);
302	}
303	}
304	}
305	else if (!strcmp(option,"All")) {
306	for (Int_t i = 0; i< numpart; i++) {
307	new(Particles[i]) TParticle(
308	fHepevt->IDHEP[i],
309	fHepevt->ISTHEP[i],
310	fHepevt->JMOHEP[i][0]-1,
311	fHepevt->JMOHEP[i][1]-1,
312	fHepevt->JDAHEP[i][0]-1,
313	fHepevt->JDAHEP[i][1]-1,
314
315	fHepevt->PHEP[i][0],
316	fHepevt->PHEP[i][1],
317	fHepevt->PHEP[i][2],
318	fHepevt->PHEP[i][3],
319	fHepevt->VHEP[i][0],
320	fHepevt->VHEP[i][1],
321	fHepevt->VHEP[i][2],
322	fHepevt->VHEP[i][3]);
323	}
324	}
325	return numpart;
326	}
327
328	void THerwig6::Hwigin()
329	{
330	hwigin_();
331	}
332
333	void THerwig6::Hwuinc()
334	{
335	hwuinc_();
336	}
337
338	void THerwig6::Hwusta(char* name)
339
340	{
341	hwusta_(name,8);
342	}
343
344	void THerwig6::Hweini()
345
346	{
347	hweini_();
348	}
349
350	void THerwig6::Hwuine()
351
352	{
353	hwuine_();
354	}
355
356	void THerwig6::Hwepro()
357
358	{
359	hwepro_();
360	}
361
362	void THerwig6::Hwbgen()
363
364	{
365	hwbgen_();
366	}
367
368	void THerwig6::Hwdhob()
369
370	{
371	hwdhob_();
372	}
373
374	void THerwig6::Hwcfor()
375
376	{
377	hwcfor_();
378	}
379
380	void THerwig6::Hwcdec()
381
382	{
383	hwcdec_();
384	}
385
386	void THerwig6::Hwdhad()
387
388	{
389	hwdhad_();
390	}
391
392	void THerwig6::Hwdhvy()
393
394	{
395	hwdhvy_();
396	}
397
398	void THerwig6::Hwmevt()
399
400	{
401	hwmevt_();
402	}
403
404	void THerwig6::Hwufne()
405
406	{
407	hwufne_();
408	}
409
410	void THerwig6::Hwefin()
411
412	{
413	hwefin_();
414	}
415
416	void THerwig6::SetupTest()
417	{
418	// exampe of running herwig and generating one event
419	// after changing some options
420	Initialize("P","PBAR",900.,900.,1500);
421	// here you can set some parameters
422	SetPTMIN(15.); // Min pt in hadronic jet production
423	SetYJMIN(-4.); // Min jet rapidity
424	SetYJMAX(4.); // Max jet rapidity
425	// after you set your wished parameters
426	// herwig can do its work
427	PrepareRun();
428	int n_ev_to_generate=1;
429	for (int i=0;i<n_ev_to_generate;i++)
430	{
431	GenerateEvent();
432	// do your stuff. For ex:
433	int n_of_par=GetNumberOfParticles(); // from TGenerator
434	for (int j=0; j<n_of_par; j++)
435	{
436	TParticle* p=GetParticle(j);
437	// here you do whatever you want with the particle
438	p->Print();
439	};
440	};
441	}
442
443
444
445
446