fad2dce669e438acd354fd41c97ada67f60d999f
[u/mrichter/AliRoot.git] / 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 \n",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(const 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