THerwig/THerwig6.cxx

   1 // definition of c++ Class THerwig6 to be used in ROOT
   2 // this is a c++ interface to the F77 Herwig6 program
   3 // author: j. g. contreras jgcn@moni.mda.cinvestav.mx
   4 // date: december 22, 2000
   5 //  Class THerwig6 is an interface to the Herwig program
   6 //
   7 // C-----------------------------------------------------------------------
   8 // C                           H E R W I G
   9 // C
  10 // C            a Monte Carlo event generator for simulating
  11 // C        +---------------------------------------------------+
  12 // C        | Hadron Emission Reactions With Interfering Gluons |
  13 // C        +---------------------------------------------------+
  14 // C I.G. Knowles(*), G. Marchesini(+), M.H. Seymour($) and B.R. Webber(#)
  15 // C-----------------------------------------------------------------------
  16 // C with Minimal Supersymmetric Standard Model Matrix Elements by
  17 // C                  S. Moretti($) and K. Odagiri($)
  18 // C-----------------------------------------------------------------------
  19 // C R parity violating Supersymmetric Decays and Matrix Elements by
  20 // C                          P. Richardson(&)
  21 // C-----------------------------------------------------------------------
  22 // C matrix element corrections to top decay and Drell-Yan type processes
  23 // C                         by G. Corcella(+)
  24 // C-----------------------------------------------------------------------
  25 // C Deep Inelastic Scattering and Heavy Flavour Electroproduction by
  26 // C                  G. Abbiendi(@) and L. Stanco(%)
  27 // C-----------------------------------------------------------------------
  28 // C and Jet Photoproduction in Lepton-Hadron Collisions by J. Chyla(~)
  29 // C-----------------------------------------------------------------------
  30 // C(*)  Department of Physics & Astronomy, University of Edinburgh
  31 // C(+)  Dipartimento di Fisica, Universita di Milano
  32 // C($)  Rutherford Appleton Laboratory
  33 // C(#)  Cavendish Laboratory, Cambridge
  34 // C(&)  Department of Physics, University of Oxford
  35 // C(@)  Dipartimento di Fisica, Universita di Bologna
  36 // C(%)  Dipartimento di Fisica, Universita di Padova
  37 // C(~)  Institute of Physics, Prague
  38 // C-----------------------------------------------------------------------
  39 // C                  Version 6.100 - 16th December 1999
  40 // C-----------------------------------------------------------------------
  41 // C Main reference:
  42 // C    G.Marchesini,  B.R.Webber,  G.Abbiendi,  I.G.Knowles,  M.H.Seymour,
  43 // C    and L.Stanco, Computer Physics Communications 67 (1992) 465.
  44 // C-----------------------------------------------------------------------
  45 // C Please send e-mail about  this program  to one of the  authors at the
  46 // C following Internet addresses:
  47 // C    I.Knowles@ed.ac.uk        Giuseppe.Marchesini@mi.infn.it
  48 // C    M.Seymour@rl.ac.uk        webber@hep.phy.cam.ac.uk
  49 // C-----------------------------------------------------------------------
  50
  51
  52 #include "THerwig6.h"
  53 #include "TClonesArray.h"
  54 #include "TParticle.h"
  55 #include "TObjArray.h"
  56
  57
  58 ClassImp(THerwig6)
  59
  60 extern "C" {
  61   void*  herwig6_common_block_address_(char*, int len);
  62   void   herwig6_open_fortran_file_ (int* lun, char* name, int);
  63   void   herwig6_close_fortran_file_(int* lun);
  64 }
  65
  66
  67 THerwig6::THerwig6() : TGenerator("Herwig6","Herwig6"),
  68   fHepevt((Hepevt_t*) herwig6_common_block_address_((char*)"HEPEVT",6)),
  69   fHwbeam((Hwbeam_t*) herwig6_common_block_address_((char*)"HWBEAM",6)),
  70   fHwbmch((Hwbmch_t*) herwig6_common_block_address_((char*)"HWBMCH",6)),
  71   fHwproc((Hwproc_t*) herwig6_common_block_address_((char*)"HWPROC",6)),
  72   fHwpram((Hwpram_t*) herwig6_common_block_address_((char*)"HWPRAM",6)),
  73   fHwprch((Hwprch_t*) herwig6_common_block_address_((char*)"HWPRCH",6)),
  74   fHwpart((Hwpart_t*) herwig6_common_block_address_((char*)"HWPART",6)),
  75   fHwparp((Hwparp_t*) herwig6_common_block_address_((char*)"HWPARP",6)),
  76   fHwbosc((Hwbosc_t*) herwig6_common_block_address_((char*)"HWBOSC",6)),
  77   fHwparc((Hwparc_t*) herwig6_common_block_address_((char*)"HWPARC",6)),
  78   fHwbrch((Hwbrch_t*) herwig6_common_block_address_((char*)"HWBRCH",6)),
  79   fHwevnt((Hwevnt_t*) herwig6_common_block_address_((char*)"HWEVNT",6)),
  80   fHwhard((Hwhard_t*) herwig6_common_block_address_((char*)"HWHARD",6)),
  81   fHwprop((Hwprop_t*) herwig6_common_block_address_((char*)"HWPROP",6)),
  82   fHwunam((Hwunam_t*) herwig6_common_block_address_((char*)"HWUNAM",6)),
  83   fHwupdt((Hwupdt_t*) herwig6_common_block_address_((char*)"HWUPDT",6)),
  84   fHwuwts((Hwuwts_t*) herwig6_common_block_address_((char*)"HWUWTS",6)),
  85   fHwuclu((Hwuclu_t*) herwig6_common_block_address_((char*)"HWUCLU",6)),
  86   fHwdist((Hwdist_t*) herwig6_common_block_address_((char*)"HWDIST",6)),
  87   fHwqdks((Hwqdks_t*) herwig6_common_block_address_((char*)"HWQDKS",6)),
  88   fHwusud((Hwusud_t*) herwig6_common_block_address_((char*)"HWUSUD",6)),
  89   fHwsusy((Hwsusy_t*) herwig6_common_block_address_((char*)"HWSUSY",6)),
  90   fHwrpar((Hwrpar_t*) herwig6_common_block_address_((char*)"HWRPAR",6)),
  91   fHwminb((Hwminb_t*) herwig6_common_block_address_((char*)"HWMINB",6)),
  92   fHwclus((Hwclus_t*) herwig6_common_block_address_((char*)"HWCLUS",6)),
  93   fHwgrav((Hwgrav_t*) herwig6_common_block_address_((char*)"HWGRAV",6)),
  94   fHw6202((Hw6202_t*) herwig6_common_block_address_((char*)"HW6202",6)),
  95   fHw6203((Hw6203_t*) herwig6_common_block_address_((char*)"HW6203",6)),
  96   fHw6300((Hw6300_t*) herwig6_common_block_address_((char*)"HW6300",6)),
  97   fHwpmrs((Hwpmrs_t*) herwig6_common_block_address_((char*)"HWPMRS",6)),
  98   fHwcirc((Hwcirc_t*) herwig6_common_block_address_((char*)"HWCIRC",6)),
  99   fHwdspb((Hwdspb_t*) herwig6_common_block_address_((char*)"HWDSPB",6)),
 100   fHwdsp2((Hwdsp2_t*) herwig6_common_block_address_((char*)"HWDSP2",6)),
 101   fHwdsp3((Hwdsp3_t*) herwig6_common_block_address_((char*)"HWDSP3",6)),
 102   fHwdsp4((Hwdsp4_t*) herwig6_common_block_address_((char*)"HWDSP4",6)),
 103   fHwdspn((Hwdspn_t*) herwig6_common_block_address_((char*)"HWDSPN",6)),
 104   fHwspin((Hwspin_t*) herwig6_common_block_address_((char*)"HWSPIN",6)),
 105   fHwstau((Hwstau_t*) herwig6_common_block_address_((char*)"HWSTAU",6)),
 106   fHwgupr((Hwgupr_t*) herwig6_common_block_address_((char*)"HWGUPR",6)),
 107   fHw6500((Hw6500_t*) herwig6_common_block_address_((char*)"HW6500",6)),
 108   fHw6504((Hw6504_t*) herwig6_common_block_address_((char*)"HW6504",6)),
 109   fHw6506((Hw6506_t*) herwig6_common_block_address_((char*)"HW6506",6))
 110 {
 111
 112 // THerwig6 constructor: creates a TClonesArray in which it will store all
 113 // particles. Note that there may be only one functional THerwig6 object
 114 // at a time, so it's not use to create more than one instance of it.
 115
 116   delete fParticles; // was allocated as TObjArray in TGenerator
 117   fParticles = new TClonesArray("TParticle",50);
 118
 119   // initialize common-blocks
 120  }
 121
 122 THerwig6::THerwig6(const THerwig6 & source): TGenerator(source),
 123   fHepevt((Hepevt_t*) herwig6_common_block_address_((char*)"HEPEVT",6)),
 124   fHwbeam((Hwbeam_t*) herwig6_common_block_address_((char*)"HWBEAM",6)),
 125   fHwbmch((Hwbmch_t*) herwig6_common_block_address_((char*)"HWBMCH",6)),
 126   fHwproc((Hwproc_t*) herwig6_common_block_address_((char*)"HWPROC",6)),
 127   fHwpram((Hwpram_t*) herwig6_common_block_address_((char*)"HWPRAM",6)),
 128   fHwprch((Hwprch_t*) herwig6_common_block_address_((char*)"HWPRCH",6)),
 129   fHwpart((Hwpart_t*) herwig6_common_block_address_((char*)"HWPART",6)),
 130   fHwparp((Hwparp_t*) herwig6_common_block_address_((char*)"HWPARP",6)),
 131   fHwbosc((Hwbosc_t*) herwig6_common_block_address_((char*)"HWBOSC",6)),
 132   fHwparc((Hwparc_t*) herwig6_common_block_address_((char*)"HWPARC",6)),
 133   fHwbrch((Hwbrch_t*) herwig6_common_block_address_((char*)"HWBRCH",6)),
 134   fHwevnt((Hwevnt_t*) herwig6_common_block_address_((char*)"HWEVNT",6)),
 135   fHwhard((Hwhard_t*) herwig6_common_block_address_((char*)"HWHARD",6)),
 136   fHwprop((Hwprop_t*) herwig6_common_block_address_((char*)"HWPROP",6)),
 137   fHwunam((Hwunam_t*) herwig6_common_block_address_((char*)"HWUNAM",6)),
 138   fHwupdt((Hwupdt_t*) herwig6_common_block_address_((char*)"HWUPDT",6)),
 139   fHwuwts((Hwuwts_t*) herwig6_common_block_address_((char*)"HWUWTS",6)),
 140   fHwuclu((Hwuclu_t*) herwig6_common_block_address_((char*)"HWUCLU",6)),
 141   fHwdist((Hwdist_t*) herwig6_common_block_address_((char*)"HWDIST",6)),
 142   fHwqdks((Hwqdks_t*) herwig6_common_block_address_((char*)"HWQDKS",6)),
 143   fHwusud((Hwusud_t*) herwig6_common_block_address_((char*)"HWUSUD",6)),
 144   fHwsusy((Hwsusy_t*) herwig6_common_block_address_((char*)"HWSUSY",6)),
 145   fHwrpar((Hwrpar_t*) herwig6_common_block_address_((char*)"HWRPAR",6)),
 146   fHwminb((Hwminb_t*) herwig6_common_block_address_((char*)"HWMINB",6)),
 147   fHwclus((Hwclus_t*) herwig6_common_block_address_((char*)"HWCLUS",6)),
 148   fHwgrav((Hwgrav_t*) herwig6_common_block_address_((char*)"HWGRAV",6)),
 149   fHw6202((Hw6202_t*) herwig6_common_block_address_((char*)"HW6202",6)),
 150   fHw6203((Hw6203_t*) herwig6_common_block_address_((char*)"HW6203",6)),
 151   fHw6300((Hw6300_t*) herwig6_common_block_address_((char*)"HW6300",6)),
 152   fHwpmrs((Hwpmrs_t*) herwig6_common_block_address_((char*)"HWPMRS",6)),
 153   fHwcirc((Hwcirc_t*) herwig6_common_block_address_((char*)"HWCIRC",6)),
 154   fHwdspb((Hwdspb_t*) herwig6_common_block_address_((char*)"HWDSPB",6)),
 155   fHwdsp2((Hwdsp2_t*) herwig6_common_block_address_((char*)"HWDSP2",6)),
 156   fHwdsp3((Hwdsp3_t*) herwig6_common_block_address_((char*)"HWDSP3",6)),
 157   fHwdsp4((Hwdsp4_t*) herwig6_common_block_address_((char*)"HWDSP4",6)),
 158   fHwdspn((Hwdspn_t*) herwig6_common_block_address_((char*)"HWDSPN",6)),
 159   fHwspin((Hwspin_t*) herwig6_common_block_address_((char*)"HWSPIN",6)),
 160   fHwstau((Hwstau_t*) herwig6_common_block_address_((char*)"HWSTAU",6)),
 161   fHwgupr((Hwgupr_t*) herwig6_common_block_address_((char*)"HWGUPR",6)),
 162   fHw6500((Hw6500_t*) herwig6_common_block_address_((char*)"HW6500",6)),
 163   fHw6504((Hw6504_t*) herwig6_common_block_address_((char*)"HW6504",6)),
 164   fHw6506((Hw6506_t*) herwig6_common_block_address_((char*)"HW6506",6))
 165 {
 166     Fatal("THerwig6","Copy constructor not implemented yet");
 167 }
 168 //------------------------------------------------------------------------------
 169  THerwig6::~THerwig6()
 170  {
 171    // Destructor. The data members of TGenerator are delete by itself
 172  }
 173
 174 //______________________________________________________________________________
 175 void THerwig6::GenerateEvent()
 176 {
 177
 178   //initialize event
 179   hwuine_();
 180   // generate hard subprocess
 181   hwepro_();
 182   // generate parton cascades
 183   hwbgen_();
 184   // do heavy objects decay
 185   hwdhob_();
 186   // do cluster formation
 187   hwcfor_();
 188   // do cluster decays
 189   hwcdec_();
 190   // do unstable particle decays
 191   hwdhad_();
 192   // do heavy flavor hadrons decay
 193   hwdhvy_();
 194   // add soft underlying event
 195   hwmevt_();
 196   // finish event
 197   hwufne_();
 198 }
 199 //______________________________________________________________________________
 200 void THerwig6::OpenFortranFile(int lun, char* name) {
 201   herwig6_open_fortran_file_(&lun, name, strlen(name));
 202 }
 203
 204 //______________________________________________________________________________
 205 void THerwig6::CloseFortranFile(int lun) {
 206   herwig6_close_fortran_file_(&lun);
 207 }
 208
 209 void THerwig6::Initialize(const char *beam, const char *target, double pbeam1, double pbeam2, int iproc)
 210
 211 {
 212   // perform the initialization for Herwig6
 213   // sets correct title.
 214   // after calling this method all parameters are set to their default
 215   // values. If you want to modify any parameter you have to set the new
 216   // value after calling Initialize and before PrepareRun.
 217
 218    char  cbeam[8];
 219    strncpy(cbeam,beam,8);
 220    char  ctarget[8];
 221    strncpy(ctarget,target,8);
 222    printf("\n Initializing Herwig !! \n");
 223    if ( (!strncmp(beam, "E+"    ,2)) &&
 224         (!strncmp(beam, "E-"    ,2)) &&
 225         (!strncmp(beam, "MU+"   ,3)) &&
 226         (!strncmp(beam, "MU-"   ,3)) &&
 227         (!strncmp(beam, "NUE"   ,3)) &&
 228         (!strncmp(beam, "NUEB"  ,4)) &&
 229         (!strncmp(beam, "NUMU"  ,4)) &&
 230         (!strncmp(beam, "NMUB"  ,4)) &&
 231         (!strncmp(beam, "NTAU"  ,4)) &&
 232         (!strncmp(beam, "NTAB"  ,4)) &&
 233         (!strncmp(beam, "GAMA"  ,4)) &&
 234         (!strncmp(beam, "P       ",8)) &&
 235         (!strncmp(beam, "PBAR    ",8)) &&
 236         (!strncmp(beam, "N"     ,1)) &&
 237         (!strncmp(beam, "NBAR"  ,4)) &&
 238         (!strncmp(beam, "PI+"   ,3)) &&
 239         (!strncmp(beam, "PI-"   ,3)) ) {
 240       printf("WARNING! In THerwig6:Initialize():\n");
 241       printf(" specified beam=%s is unrecognized .\n",beam);
 242       printf(" resetting to \"P\" .");
 243       sprintf(cbeam,"P");
 244    }
 245
 246    if ( (!strncmp(target, "E+"    ,2)) &&
 247         (!strncmp(target, "E-"    ,2)) &&
 248         (!strncmp(target, "MU+"   ,3)) &&
 249         (!strncmp(target, "MU-"   ,3)) &&
 250         (!strncmp(target, "NUE"   ,3)) &&
 251         (!strncmp(target, "NUEB"  ,4)) &&
 252         (!strncmp(target, "NUMU"  ,4)) &&
 253         (!strncmp(target, "NMUB"  ,4)) &&
 254         (!strncmp(target, "NTAU"  ,4)) &&
 255         (!strncmp(target, "NTAB"  ,4)) &&
 256         (!strncmp(target, "GAMA"  ,4)) &&
 257         (!strncmp(target, "P       ",8)) &&
 258         (!strncmp(target, "PBAR    ",8)) &&
 259         (!strncmp(target, "N"     ,1)) &&
 260         (!strncmp(target, "NBAR"  ,4)) &&
 261         (!strncmp(target, "PI+"   ,3)) &&
 262         (!strncmp(target, "PI-"   ,3)) ) {
 263       printf("WARNING! In THerwig6:Initialize():\n");
 264       printf(" specified target=%s is unrecognized .\n",target);
 265       printf(" resetting to \"P\" .");
 266       sprintf(ctarget,"P");
 267    }
 268
 269    // initialization:
 270    // type of beams
 271    strncpy(fHwbmch->PART1,beam,8);
 272    strncpy(fHwbmch->PART2,target,8);
 273    // momentum of beams
 274    fHwproc->PBEAM1=pbeam1;
 275    fHwproc->PBEAM2=pbeam2;
 276    // process to generate
 277    fHwproc->IPROC=iproc;
 278    // not used in the class definition
 279    fHwproc->MAXEV=1;
 280
 281    // reset all parameters
 282    hwigin_();
 283
 284    // set correct title
 285    char atitle[132];
 286    double win=pbeam1+pbeam2;
 287    printf("\n %s - %s at %g GeV",beam,target,win);
 288    sprintf(atitle,"%s-%s at %g GeV",cbeam,ctarget,win);
 289    SetTitle(atitle);
 290 }
 291
 292 void THerwig6::InitializeJimmy(const char *beam, const char *target, double pbeam1, double pbeam2, int iproc)
 293
 294 {
 295   // perform the initialization for Herwig6
 296   // sets correct title.
 297   // after calling this method all parameters are set to their default
 298   // values. If you want to modify any parameter you have to set the new
 299   // value after calling Initialize and before PrepareRun.
 300
 301    char  cbeam[8];
 302    strncpy(cbeam,beam,8);
 303    char  ctarget[8];
 304    strncpy(ctarget,target,8);
 305    printf("\n Initializing Herwig !! \n");
 306    if ( (!strncmp(beam, "E+"    ,2)) &&
 307         (!strncmp(beam, "E-"    ,2)) &&
 308         (!strncmp(beam, "MU+"   ,3)) &&
 309         (!strncmp(beam, "MU-"   ,3)) &&
 310         (!strncmp(beam, "NUE"   ,3)) &&
 311         (!strncmp(beam, "NUEB"  ,4)) &&
 312         (!strncmp(beam, "NUMU"  ,4)) &&
 313         (!strncmp(beam, "NMUB"  ,4)) &&
 314         (!strncmp(beam, "NTAU"  ,4)) &&
 315         (!strncmp(beam, "NTAB"  ,4)) &&
 316         (!strncmp(beam, "GAMA"  ,4)) &&
 317         (!strncmp(beam, "P       ",8)) &&
 318         (!strncmp(beam, "PBAR    ",8)) &&
 319         (!strncmp(beam, "N"     ,1)) &&
 320         (!strncmp(beam, "NBAR"  ,4)) &&
 321         (!strncmp(beam, "PI+"   ,3)) &&
 322         (!strncmp(beam, "PI-"   ,3)) ) {
 323       printf("WARNING! In THerwig6:Initialize():\n");
 324       printf(" specified beam=%s is unrecognized .\n",beam);
 325       printf(" resetting to \"P\" .");
 326       sprintf(cbeam,"P");
 327    }
 328
 329    if ( (!strncmp(target, "E+"    ,2)) &&
 330         (!strncmp(target, "E-"    ,2)) &&
 331         (!strncmp(target, "MU+"   ,3)) &&
 332         (!strncmp(target, "MU-"   ,3)) &&
 333         (!strncmp(target, "NUE"   ,3)) &&
 334         (!strncmp(target, "NUEB"  ,4)) &&
 335         (!strncmp(target, "NUMU"  ,4)) &&
 336         (!strncmp(target, "NMUB"  ,4)) &&
 337         (!strncmp(target, "NTAU"  ,4)) &&
 338         (!strncmp(target, "NTAB"  ,4)) &&
 339         (!strncmp(target, "GAMA"  ,4)) &&
 340         (!strncmp(target, "P       ",8)) &&
 341         (!strncmp(target, "PBAR    ",8)) &&
 342         (!strncmp(target, "N"     ,1)) &&
 343         (!strncmp(target, "NBAR"  ,4)) &&
 344         (!strncmp(target, "PI+"   ,3)) &&
 345         (!strncmp(target, "PI-"   ,3)) ) {
 346       printf("WARNING! In THerwig6:Initialize():\n");
 347       printf(" specified target=%s is unrecognized .\n",target);
 348       printf(" resetting to \"P\" .");
 349       sprintf(ctarget,"P");
 350    }
 351
 352    // initialization:
 353    // type of beams
 354    strncpy(fHwbmch->PART1,beam,8);
 355    strncpy(fHwbmch->PART2,target,8);
 356    // momentum of beams
 357    fHwproc->PBEAM1=pbeam1;
 358    fHwproc->PBEAM2=pbeam2;
 359    // process to generate
 360    fHwproc->IPROC=iproc;
 361    // not used in the class definition
 362    fHwproc->MAXEV=1;
 363
 364    // reset all parameters
 365    hwigin_();
 366    // JIMMY initialization
 367    jimmin_();
 368
 369    // set correct title
 370    char atitle[132];
 371    double win=pbeam1+pbeam2;
 372    printf("\n %s - %s at %g GeV",beam,target,win);
 373    sprintf(atitle,"%s-%s at %g GeV",cbeam,ctarget,win);
 374    SetTitle(atitle);
 375 }
 376
 377 void THerwig6::PrepareRun()
 378 {
 379   // compute parameter dependent constants
 380   hwuinc_();
 381   // initialize elementary processes
 382   hweini_();
 383 }
 384
 385 void THerwig6::PrepareRunJimmy()
 386 {
 387   // compute parameter dependent constants
 388   hwuinc_();
 389   // initialize elementary processes
 390   hweini_();
 391   // more initializations for JIMMY
 392   jminit_();
 393 }
 394 //______________________________________________________________________________
 395 TObjArray* THerwig6::ImportParticles(Option_t *option)
 396 {
 397 //
 398 //  Default primary creation method. It reads the /HEPEVT/ common block which
 399 //  has been filled by the GenerateEvent method. If the event generator does
 400 //  not use the HEPEVT common block, This routine has to be overloaded by
 401 //  the subclasses.
 402 //  The default action is to store only the stable particles (ISTHEP = 1)
 403 //  This can be demanded explicitly by setting the option = "Final"
 404 //  If the option = "All", all the particles are stored.
 405 //
 406   fParticles->Clear();
 407   Int_t numpart = fHepevt->NHEP;
 408   TClonesArray &a = *((TClonesArray*)fParticles);
 409   if (!strcmp(option,"") || !strcmp(option,"Final")) {
 410     for (Int_t i = 0; i<=numpart; i++) {
 411       if (fHepevt->ISTHEP[i] == 1) {
 412 //
 413 //  Use the common block values for the TParticle constructor
 414 //
 415         new(a[i]) TParticle(
 416                                    fHepevt->IDHEP[i],
 417                                    fHepevt->ISTHEP[i],
 418                                    fHepevt->JMOHEP[i][0]-1,
 419                                    fHepevt->JMOHEP[i][1]-1,
 420                                    fHepevt->JDAHEP[i][0]-1,
 421                                    fHepevt->JDAHEP[i][1]-1,
 422
 423                                    fHepevt->PHEP[i][0],
 424                                    fHepevt->PHEP[i][1],
 425                                    fHepevt->PHEP[i][2],
 426                                    fHepevt->PHEP[i][3],
 427                                    fHepevt->VHEP[i][0],
 428                                    fHepevt->VHEP[i][1],
 429                                    fHepevt->VHEP[i][2],
 430                                    fHepevt->VHEP[i][3]);
 431         }
 432      }
 433   }
 434   else if (!strcmp(option,"All")) {
 435     for (Int_t i = 0; i<=numpart; i++) {
 436       new(a[i]) TParticle(
 437                                    fHepevt->IDHEP[i],
 438                                    fHepevt->ISTHEP[i],
 439                                    fHepevt->JMOHEP[i][0]-1,
 440                                    fHepevt->JMOHEP[i][1]-1,
 441                                    fHepevt->JDAHEP[i][0]-1,
 442                                    fHepevt->JDAHEP[i][1]-1,
 443
 444                                    fHepevt->PHEP[i][0],
 445                                    fHepevt->PHEP[i][1],
 446                                    fHepevt->PHEP[i][2],
 447                                    fHepevt->PHEP[i][3],
 448                                    fHepevt->VHEP[i][0],
 449                                    fHepevt->VHEP[i][1],
 450                                    fHepevt->VHEP[i][2],
 451                                    fHepevt->VHEP[i][3]);
 452     }
 453   }
 454   return fParticles;
 455 }
 456
 457 //______________________________________________________________________________
 458 Int_t THerwig6::ImportParticles(TClonesArray *particles, Option_t *option)
 459 {
 460 //
 461 //  Default primary creation method. It reads the /HEPEVT/ common block which
 462 //  has been filled by the GenerateEvent method. If the event generator does
 463 //  not use the HEPEVT common block, This routine has to be overloaded by
 464 //  the subclasses.
 465 //  The function loops on the generated particles and store them in
 466 //  the TClonesArray pointed by the argument particles.
 467 //  The default action is to store only the stable particles (ISTHEP = 1)
 468 //  This can be demanded explicitly by setting the option = "Final"
 469 //  If the option = "All", all the particles are stored.
 470 //
 471   if (particles == 0) return 0;
 472   TClonesArray &refParticles = *particles;
 473   refParticles.Clear();
 474   Int_t numpart = fHepevt->NHEP;
 475   if (!strcmp(option,"") || !strcmp(option,"Final")) {
 476     for (Int_t i = 0; i< numpart; i++) {
 477       if (fHepevt->ISTHEP[i] == 1) {
 478 //
 479 //  Use the common block values for the TParticle constructor
 480 //
 481         new(refParticles[i]) TParticle(
 482                                    fHepevt->IDHEP[i],
 483                                    fHepevt->ISTHEP[i],
 484                                    fHepevt->JMOHEP[i][0]-1,
 485                                    fHepevt->JMOHEP[i][1]-1,
 486                                    fHepevt->JDAHEP[i][0]-1,
 487                                    fHepevt->JDAHEP[i][1]-1,
 488
 489                                    fHepevt->PHEP[i][0],
 490                                    fHepevt->PHEP[i][1],
 491                                    fHepevt->PHEP[i][2],
 492                                    fHepevt->PHEP[i][3],
 493                                    fHepevt->VHEP[i][0],
 494                                    fHepevt->VHEP[i][1],
 495                                    fHepevt->VHEP[i][2],
 496                                    fHepevt->VHEP[i][3]);
 497         }
 498      }
 499   }
 500   else if (!strcmp(option,"All")) {
 501     for (Int_t i = 0; i< numpart; i++) {
 502       new(refParticles[i]) TParticle(
 503                                    fHepevt->IDHEP[i],
 504                                    fHepevt->ISTHEP[i],
 505                                    fHepevt->JMOHEP[i][0]-1,
 506                                    fHepevt->JMOHEP[i][1]-1,
 507                                    fHepevt->JDAHEP[i][0]-1,
 508                                    fHepevt->JDAHEP[i][1]-1,
 509
 510                                    fHepevt->PHEP[i][0],
 511                                    fHepevt->PHEP[i][1],
 512                                    fHepevt->PHEP[i][2],
 513                                    fHepevt->PHEP[i][3],
 514                                    fHepevt->VHEP[i][0],
 515                                    fHepevt->VHEP[i][1],
 516                                    fHepevt->VHEP[i][2],
 517                                    fHepevt->VHEP[i][3]);
 518     }
 519   }
 520   return numpart;
 521 }
 522
 523 void THerwig6::Hwigin()
 524 {
 525   hwigin_();
 526 }
 527
 528 void THerwig6::Hwuinc()
 529 {
 530   hwuinc_();
 531 }
 532
 533 void THerwig6::Hwusta(char* name)
 534
 535 {
 536   hwusta_(name,8);
 537 }
 538
 539 void THerwig6::Hweini()
 540
 541 {
 542   hweini_();
 543 }
 544
 545 void THerwig6::Hwuine()
 546
 547 {
 548   hwuine_();
 549 }
 550
 551 void THerwig6::Hwepro()
 552
 553 {
 554   hwepro_();
 555 }
 556
 557 void THerwig6::Hwbgen()
 558
 559 {
 560   hwbgen_();
 561 }
 562
 563 void THerwig6::Hwdhob()
 564
 565 {
 566   hwdhob_();
 567 }
 568
 569 void THerwig6::Hwcfor()
 570
 571 {
 572   hwcfor_();
 573 }
 574
 575 void THerwig6::Hwcdec()
 576
 577 {
 578   hwcdec_();
 579 }
 580
 581 void THerwig6::Hwdhad()
 582
 583 {
 584   hwdhad_();
 585 }
 586
 587 void THerwig6::Hwdhvy()
 588
 589 {
 590   hwdhvy_();
 591 }
 592
 593 void THerwig6::Hwmevt()
 594
 595 {
 596   hwmevt_();
 597 }
 598
 599 void THerwig6::Hwufne()
 600
 601 {
 602   hwufne_();
 603 }
 604
 605 void THerwig6::Hwefin()
 606
 607 {
 608   hwefin_();
 609 }
 610
 611 void THerwig6::Hwiodk(int iopt)
 612
 613 {
 614   hwiodk_(iopt);
 615 }
 616
 617 void THerwig6::SetupTest()
 618 {
 619   // exampe of running herwig and generating one event
 620   // after changing some options
 621   Initialize("P","PBAR",900.,900.,1500);
 622   // here you can set some parameters
 623   SetPTMIN(15.); // Min pt in hadronic jet production
 624   SetYJMIN(-4.); // Min jet rapidity
 625   SetYJMAX(4.);  // Max jet rapidity
 626   // after you set your wished parameters
 627   // herwig can do its work
 628   PrepareRun();
 629   int nEvToGenerate=1;
 630   for (int i=0;i<nEvToGenerate;i++)
 631     {
 632       GenerateEvent();
 633       // do your stuff. For ex:
 634       int nOfPar=GetNumberOfParticles(); // from TGenerator
 635       for (int j=0; j<nOfPar; j++)
 636         {
 637           TParticle* p=GetParticle(j);
 638           // here you do whatever you want with the particle
 639           p->Print();
 640         };
 641     };
 642 }
 643
 644 // Jimmy subroutines
 645
 646 void THerwig6::Jminit()
 647 {
 648   jminit_();
 649 }
 650
 651 void THerwig6::Jimmin()
 652 {
 653   jimmin_();
 654 }
 655
 656 void THerwig6::Jmefin()
 657 {
 658   jmefin_();
 659 }
 660
 661
 662
 663
 664
 665