Bringing CMakeLists under svn maintenance
[u/mrichter/AliRoot.git] / THijing / THijing.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /* $Id$ */
17
18 ////////////////////////////////////////////////////////////////////////////////
19 //                                                                            //
20 // THijing                                                                    //
21 //                                                                            //
22 // THijing is an interface class to F77 version of Hijing 1.36                //
23 // event generator, written by X.N. Wang and  M. Gyulassy.                    //
24 // For details see http://nta2.lbl.gov/~xnwang                                //
25 //                                                                            //
26 //          **************************************************                //
27 //          *     |             _______      /  ------/      *                //
28 //          *   ----- ------     |_____|     /_/     /       *                //
29 //          *    ||    /        |_____|      /    /          *                //
30 //          *    /|   /_/       /_______    /_  /    _       *                //
31 //          *   / |     / /     /  /  / |        -------     *                //
32 //          *     |    / /       /  /  |     /     |         *                //
33 //          *     |   / /       /  / _|    /   -------       *                //
34 //          *                                                *                //
35 //          **************************************************                //
36 //                                HIJING                                      //
37 //                 Heavy Ion Jet INteraction Generator                        //
38 //                                  by                                        //
39 //                     X. N. Wang  and  M. Gyulassy                           //
40 //                     Lawrence Berkeley Laboratory                           //
41 //****************************************************************************//
42
43
44 #include <TClonesArray.h>
45 #include <TObjArray.h>
46 #include <TParticle.h>
47 #include <TROOT.h>
48
49 #include "Hcommon.h"
50 #include "THijing.h"
51
52 #ifndef WIN32
53 # define hijset hijset_
54 # define hijing hijing_
55 # define profile profile_
56 # define rluget_hijing rluget_hijing_
57 # define rluset_hijing rluset_hijing_
58 # define lulist_hijing lulist_hijing_
59 # define type_of_call
60 #else
61 # define hijset HIJSET
62 # define hijing HIJING
63 # define profile PROFILE
64 # define rluget_hijing RLUGET_HIJING
65 # define rluset_hijing RLUSET_HIJING
66 # define lulist_hijing LULIST_HIJING
67 # define type_of_call _stdcall
68 #endif
69
70 #ifndef WIN32
71 extern "C" void type_of_call hijset(Float_t & , const char *, 
72                                    const char *, const char *,
73                                    Int_t & , Int_t &, Int_t &,
74                                    Int_t &, const int, 
75                                    const int, const int);
76 extern "C" float type_of_call profile(Float_t &);
77 extern "C" void type_of_call hijing(const char *, Float_t  &,
78                                    Float_t &, const int);
79 extern "C" void type_of_call rluget_hijing(Int_t & lfn, Int_t & move);
80
81 extern "C" void type_of_call rluset_hijing(Int_t & lfn, Int_t & move);
82 extern "C" void type_of_call lulist_hijing(Int_t &);
83
84 #else
85 #endif
86
87
88
89 ClassImp(THijing)
90
91
92 THijing::THijing(): 
93     TGenerator("Hijing","Hijing"),
94     fEfrm(5500.),
95     fFrame("CMS"),
96     fProj("A"),
97     fTarg("A"),
98     fIap(208),
99     fIzp(82),
100     fIat(208),
101     fIzt(82),
102     fBmin(0.),
103     fBmax(5.)
104 {
105 // Default constructor 
106 }
107
108 //______________________________________________________________________________
109 THijing::THijing(Float_t efrm, const char *frame="CMS", 
110                  const char *proj="A", const char *targ="A", 
111                  Int_t iap=207, Int_t izp=82, Int_t iat=207, Int_t izt=82,
112                  Float_t bmin=0, Float_t bmax=20): 
113     TGenerator("Hijing","Hijing"),
114     fEfrm(efrm),
115     fFrame(frame),
116     fProj(proj),
117     fTarg(targ),
118     fIap(iap),
119     fIzp(izp),
120     fIat(iat),
121     fIzt(izt),
122     fBmin(bmin),
123     fBmax(bmax)
124 {
125 // THijing constructor: 
126 // Note that there may be only one functional THijing object
127 // at a time, so it's not use to create more than one instance of it.
128 }
129
130 //______________________________________________________________________________
131 THijing::~THijing()
132 {
133 // Destructor
134 }
135
136
137 TObjArray* THijing::ImportParticles(Option_t *option)
138 {
139 //
140 //  Default primary creation method. It reads the /HEPEVT/ common block which
141 //  has been filled by the GenerateEvent method. If the event generator does
142 //  not use the HEPEVT common block, This routine has to be overloaded by
143 //  the subclasses.
144 //  The function loops on the generated particles and store them in
145 //  the TClonesArray pointed by the argument particles.
146 //  The default action is to store only the stable particles (ISTHEP = 1)
147 //  This can be demanded explicitly by setting the option = "Final"
148 //  If the option = "All", all the particles are stored.
149 //
150     fParticles->Clear();
151     Int_t numpart = HIMAIN1.natt;
152     printf("\n THijing: HIJING stack contains %d particles.", numpart);
153     printf("\n THijing: Total energy:         %f           ", HIMAIN1.eatt);
154     printf("\n THijing: Number of hard scatterings: %d     ", HIMAIN1.jatt);
155     Int_t nump = 0;
156     if (!strcmp(option,"") || !strcmp(option,"Final")) {
157         for (Int_t i = 0; i < numpart; i++) {
158           
159             if (HIMAIN2.katt[3][i] == 1) {
160 //
161 //  Use the common block values for the TParticle constructor
162 //
163                 nump++;
164                 TParticle* p = new TParticle(
165                     HIMAIN2.katt[0][i], HIMAIN2.katt[3][i] ,
166                     -1, -1, -1, -1,
167                     HIMAIN2.patt[0][i], HIMAIN2.patt[1][i], HIMAIN2.patt[2][i], HIMAIN2.patt[3][i] ,
168                     HIMAIN2.vatt[0][i], HIMAIN2.vatt[1][i], HIMAIN2.vatt[2][i], HIMAIN2.vatt[3][i] 
169                   );
170                 p->SetUniqueID(HIMAIN2.katt[1][i]);
171                 fParticles->Add(p);
172             }
173         }
174     }
175     else if (!strcmp(option,"All")) {
176         nump = numpart; 
177         for (Int_t i = 0; i < numpart; i++) {
178             
179             Int_t iParent = HIMAIN2.katt[2][i]-1;
180             
181             if (iParent >= 0) {
182                 TParticle *mother = (TParticle*) (fParticles->UncheckedAt(iParent));       
183                 mother->SetLastDaughter(i);
184                 if (mother->GetFirstDaughter()==-1)
185                     mother->SetFirstDaughter(i);
186             }
187             
188             TParticle* p = new TParticle(
189                 HIMAIN2.katt[0][i], HIMAIN2.katt[3][i], iParent,
190                 -1, -1, -1,
191                 HIMAIN2.patt[0][i], HIMAIN2.patt[1][i], HIMAIN2.patt[2][i], HIMAIN2.patt[3][i] ,
192                 HIMAIN2.vatt[0][i], HIMAIN2.vatt[1][i], HIMAIN2.vatt[2][i], HIMAIN2.vatt[3][i]
193                 );
194             p->SetUniqueID(HIMAIN2.katt[1][i]);
195             fParticles->Add(p);
196         }
197     }
198     return fParticles;
199 }
200
201 Int_t THijing::ImportParticles(TClonesArray *particles, Option_t *option)
202 {
203 //
204 //  Default primary creation method. It reads the /HEPEVT/ common block which
205 //  has been filled by the GenerateEvent method. If the event generator does
206 //  not use the HEPEVT common block, This routine has to be overloaded by
207 //  the subclasses.
208 //  The function loops on the generated particles and store them in
209 //  the TClonesArray pointed by the argument particles.
210 //  The default action is to store only the stable particles (ISTHEP = 1)
211 //  This can be demanded explicitly by setting the option = "Final"
212 //  If the option = "All", all the particles are stored.
213 //
214   if (particles == 0) return 0;
215   TClonesArray &particlesR = *particles;
216   particlesR.Clear();
217   Int_t numpart = HIMAIN1.natt;
218   printf("\n THijing: HIJING stack contains %d particles.", numpart);
219   printf("\n THijing: Total energy:         %f           ", HIMAIN1.eatt);
220   printf("\n THijing: Number of hard scatterings: %d     ", HIMAIN1.jatt);
221   Int_t nump = 0;
222   if (!strcmp(option,"") || !strcmp(option,"Final")) {
223       for (Int_t i = 0; i < numpart; i++) {
224           
225           if (HIMAIN2.katt[3][i] == 1) {
226 //
227 //  Use the common block values for the TParticle constructor
228 //
229             nump++;
230             new(particlesR[i]) TParticle(
231                   HIMAIN2.katt[0][i] ,
232                   HIMAIN2.katt[3][i] ,
233                   -1 ,
234                   -1,
235                   -1,
236                   -1,
237                   
238                   HIMAIN2.patt[0][i] ,
239                   HIMAIN2.patt[1][i] ,
240                   HIMAIN2.patt[2][i] ,
241                   HIMAIN2.patt[3][i] ,
242
243                   HIMAIN2.vatt[0][i] ,
244                   HIMAIN2.vatt[1][i] ,
245                   HIMAIN2.vatt[2][i] ,
246                   HIMAIN2.vatt[3][i] 
247                   );
248             particlesR[i]->SetUniqueID(HIMAIN2.katt[1][i]);
249           }
250       }
251   }
252   else if (!strcmp(option,"All")) {
253       nump = numpart; 
254       for (Int_t i = 0; i < numpart; i++) {
255
256           Int_t iParent = HIMAIN2.katt[2][i]-1;
257           
258           if (iParent >= 0) {
259               TParticle *mother = (TParticle*) (particlesR.UncheckedAt(iParent));          
260               mother->SetLastDaughter(i);
261               if (mother->GetFirstDaughter()==-1)
262                   mother->SetFirstDaughter(i);
263           }
264
265           new(particlesR[i]) TParticle(
266               HIMAIN2.katt[0][i] ,
267               HIMAIN2.katt[3][i] ,
268               iParent,
269               -1,
270               -1,
271               -1,
272               
273               HIMAIN2.patt[0][i] ,
274               HIMAIN2.patt[1][i] ,
275               HIMAIN2.patt[2][i] ,
276               HIMAIN2.patt[3][i] ,
277               
278               HIMAIN2.vatt[0][i] ,
279               HIMAIN2.vatt[1][i] ,
280               HIMAIN2.vatt[2][i] ,
281               HIMAIN2.vatt[3][i]
282               );
283           particlesR[i]->SetUniqueID(HIMAIN2.katt[1][i]);
284       }
285   }
286   return nump;
287 }
288
289 //______________________________________________________________________________
290 void THijing::SetEFRM(Float_t efrm)
291 {
292 // Set the centre of mass (CMS) or lab-energy (LAB)
293    fEfrm=efrm;
294
295 //______________________________________________________________________________
296 void THijing::SetFRAME(const char* frame)
297 {
298 // Set the frame type ("CMS" or "LAB")
299    fFrame=frame;
300
301 //______________________________________________________________________________
302 void THijing::SetPROJ(const char* proj)
303 {
304 // Set the projectile type
305    fProj=proj;
306
307 //______________________________________________________________________________
308 void THijing::SetTARG(const char* targ)
309 {
310 // Set the target type
311    fTarg=targ;
312
313 //______________________________________________________________________________
314 void THijing::SetIAP(Int_t iap)
315 {
316 // Set the projectile atomic number
317    fIap=iap;
318
319 //______________________________________________________________________________
320 void THijing::SetIZP(Int_t izp)
321 {
322 // Set the projectile charge number
323    fIzp=izp;
324
325 //______________________________________________________________________________
326 void THijing::SetIAT(Int_t iat)
327 {
328 // Set the target atomic number
329    fIat=iat;
330
331 //______________________________________________________________________________
332 void THijing::SetIZT(Int_t izt)
333 {
334 // Set the target charge number
335    fIzt=izt;
336
337 //______________________________________________________________________________
338 void THijing::SetBMIN(Float_t bmin)
339 {
340 // Set the minimum impact parameter
341    fBmin=bmin;
342
343 //______________________________________________________________________________
344 void THijing::SetBMAX(Float_t bmax)
345 {
346 // Set the maximum impact parameter
347    fBmax=bmax;
348
349 //______________________________________________________________________________
350 Float_t THijing::GetEFRM() const
351 {
352 // Get the centre of mass (CMS) or lab-energy (LAB)
353    return fEfrm;
354
355 //______________________________________________________________________________
356 const char* THijing::GetFRAME() const
357 {
358 // Get the frame type ("CMS" or "LAB")
359    return fFrame.Data();
360
361 //______________________________________________________________________________
362 const char* THijing::GetPROJ() const
363 {
364 // Get the projectile type
365    return fProj;
366
367 //______________________________________________________________________________
368 const char* THijing::GetTARG() const
369 {
370 // Set the target type
371    return fTarg;
372
373 //______________________________________________________________________________
374 Int_t THijing::GetIAP() const
375 {
376 // Get the projectile atomic number
377    return fIap;
378
379 //______________________________________________________________________________
380 Int_t THijing::GetIZP() const
381 {
382 // Get the projectile charge number
383    return fIzp;
384
385 //______________________________________________________________________________
386 Int_t THijing::GetIAT() const
387 {
388 // Get the target atomic number
389    return fIat;
390
391 //______________________________________________________________________________
392 Int_t THijing::GetIZT() const
393 {
394 // Get the target charge number
395    return fIzt;
396
397 //______________________________________________________________________________
398 Float_t THijing::GetBMIN() const
399 {
400 // Get the minimum impact parameter
401    return fBmin;
402
403 //______________________________________________________________________________
404 Float_t THijing::GetBMAX() const
405 {
406 // Get the maximum impact parameter
407    return fBmax;
408
409
410 //====================== access to common HIPARNT ===============================
411
412 //______________________________________________________________________________
413 void THijing::SetHIPR1(Int_t key,Float_t value)
414 {
415 // Set the values of array HIPR1 in common HIPARNT
416    if ( key<1 || key>100 ) {
417       printf ("ERROR in THijing:SetHIPR1(key,value): \n ");
418       printf ("      key=%i is out of range [1..100]!\n",key);
419       return;
420    }
421
422    HIPARNT.hipr1[key-1]=value;
423
424 }
425
426 //______________________________________________________________________________
427 Float_t THijing::GetHIPR1(Int_t key) const
428 {
429 // Get the values of array HIPR1 in common HIPARNT
430    if ( key<1 || key>100 ) {
431       printf ("ERROR in THijing:GetHIPR1(key): \n ");
432       printf ("      key=%i is out of range [1..100]!\n",key);
433       return 0;
434    }
435
436    return HIPARNT.hipr1[key-1];
437
438 }
439
440 //______________________________________________________________________________
441 void THijing::SetIHPR2(Int_t key,Int_t value)
442 {
443 // Set the values of array HIPR2 in common HIPARNT
444    if ( key<1 || key>50 ) {
445       printf ("ERROR in THijing:SetIHPR2(key,value): \n ");
446       printf ("      key=%i is out of range [1..50]!\n",key);
447       return;
448    }
449
450    HIPARNT.ihpr2[key-1]=value;
451
452 }
453
454 //______________________________________________________________________________
455 Int_t THijing::GetIHPR2(Int_t key) const
456 {
457 // Get the values of array HIPR2 in common HIPARNT
458    if ( key<1 || key>50 ) {
459       printf ("ERROR in THijing:GetIHPR2(key): \n ");
460       printf ("      key=%i is out of range [1..50]!\n",key);
461       return 0;
462    }
463
464    return HIPARNT.ihpr2[key-1];
465
466 }
467
468
469 //______________________________________________________________________________
470 Float_t THijing::GetHINT1(Int_t key) const
471 {
472 // Get the values of array HINT1 in common HIPARNT
473    if ( key<1 || key>100 ) {
474       printf ("ERROR in THijing:GetHINT1(key): \n ");
475       printf ("      key=%i is out of range [1..100]!\n",key);
476       return 0;
477    }
478
479    return HIPARNT.hint1[key-1];
480
481 }
482
483
484 //______________________________________________________________________________
485 Int_t THijing::GetIHNT2(Int_t key) const
486 {
487 // Get the values of array HINT2 in common HIPARNT
488    if ( key<1 || key>50 ) {
489       printf ("ERROR in THijing:GetIHNT2(key): \n ");
490       printf ("      key=%i is out of range [1..50]!\n",key);
491       return 0;
492    }
493
494    return HIPARNT.ihnt2[key-1];
495
496 }
497
498
499 //====================== access to common HIMAIN1 ===============================
500
501 //______________________________________________________________________________
502 Int_t  THijing::GetNATT() const
503 {
504 // Get the number of particles produces
505    return HIMAIN1.natt;
506
507 }
508
509 //______________________________________________________________________________
510 Float_t  THijing::GetEATT() const
511 {
512 // Get total energy of particles
513
514    return HIMAIN1.eatt;
515
516 }
517
518 //______________________________________________________________________________
519 Int_t  THijing::GetJATT() const
520 {
521 // Get number of hard scatterings
522
523    return HIMAIN1.jatt;
524
525 }
526
527 //______________________________________________________________________________
528 Int_t  THijing::GetNT() const
529 {
530 // Get number of target participants
531
532    return HIMAIN1.nt;
533
534 }
535
536 //______________________________________________________________________________
537 Int_t  THijing::GetNP() const
538 {
539 // Get number of projectile participants
540    return HIMAIN1.np;
541
542 }
543
544
545 //______________________________________________________________________________
546 Int_t  THijing::GetN0() const
547 {
548 // Get number of N-N collisions
549    return HIMAIN1.n0;
550
551 }
552 //______________________________________________________________________________
553 Int_t  THijing::GetN01() const
554 {
555 // Get number of N-wounded collisions
556
557    return HIMAIN1.n01;
558
559 }
560
561 //______________________________________________________________________________
562 Int_t  THijing::GetN10() const
563 {
564 // Get number of wounded-N collisions
565
566    return HIMAIN1.n10;
567
568 }
569
570 //______________________________________________________________________________
571 Int_t  THijing::GetN11() const
572 {
573 // Get number of wounded-wounded collisions
574
575    return HIMAIN1.n11;
576
577 }
578
579 //______________________________________________________________________________
580 Float_t  THijing::GetBB() const
581 {
582 // Get impact parameter
583
584    return HIMAIN1.bb;
585
586 }
587
588 //====================== access to common HIMAIN2 ===============================
589
590 //______________________________________________________________________________
591 Int_t THijing::GetKATT(Int_t key1, Int_t key2) const
592 {
593 // Get values of array KATT in common HIMAIN2
594    if ( key1<1 || key1>200000 ) {
595       printf("ERROR in THijing::GetKATT(key1,key2):\n");
596       printf("      key1=%i is out of range [1..200000]\n",key1);
597       return 0;
598    }
599
600    if ( key2<1 || key2>4 ) {
601       printf("ERROR in THijing::GetKATT(key1,key2):\n");
602       printf("      key2=%i is out of range [1..4]\n",key2);
603       return 0;
604    }
605    
606    return   HIMAIN2.katt[key2-1][key1-1];
607 }
608
609 //______________________________________________________________________________
610 Float_t THijing::GetPATT(Int_t key1, Int_t key2) const
611 {
612 // Get values of array PATT in common HIMAIN2
613    if ( key1<1 || key1>200000 ) {
614       printf("ERROR in THijing::GetPATT(key1,key2):\n");
615       printf("      key1=%i is out of range [1..130000]\n",key1);
616       return 0;
617    }
618
619    if ( key2<1 || key2>4 ) {
620       printf("ERROR in THijing::GetPATT(key1,key2):\n");
621       printf("      key2=%i is out of range [1..4]\n",key2);
622       return 0;
623    }
624    
625    return   HIMAIN2.patt[key2-1][key1-1];
626 }
627
628 Float_t THijing::GetVATT(Int_t key1, Int_t key2) const
629 {
630 // Get values of array VATT in common HIMAIN2
631    if ( key1<1 || key1>200000 ) {
632       printf("ERROR in THijing::GetVATT(key1,key2):\n");
633       printf("      key1=%i is out of range [1..130000]\n",key1);
634       return 0;
635    }
636
637    if ( key2<1 || key2>4 ) {
638       printf("ERROR in THijing::GetVATT(key1,key2):\n");
639       printf("      key2=%i is out of range [1..4]\n",key2);
640       return 0;
641    }
642    
643    return   HIMAIN2.vatt[key2-1][key1-1];
644 }
645
646 //====================== access to common HIJJET1 ===============================
647
648 //______________________________________________________________________________
649 Int_t THijing::GetNPJ(Int_t key) const
650 {
651 // Get values of array NPJ of common HIJJET1
652    if ( key<1 || key>300 ) {
653       printf("ERROR in THijing::GetNPJ(key):\n");
654       printf("      key=%i is out of range [1..300]\n",key);
655       return 0;
656    }
657    return HIJJET1.npj[key-1];
658 }
659
660 //______________________________________________________________________________
661 Int_t THijing::GetKFPJ(Int_t key1, Int_t key2) const
662 {
663 // Get values of array KFPJ in common HIJJET1
664    if ( key1<1 || key1>300 ) {
665       printf("ERROR in THijing::GetKFPJ(key1):\n");
666       printf("      key1=%i is out of range [1..300]\n",key1);
667       return 0;
668    }
669    if ( key2<1 || key2>500 ) {
670       printf("ERROR in THijing::GetKFPJ(key1,key2):\n");
671       printf("      key2=%i is out of range [1..500]\n",key2);
672       return 0;
673    }
674    
675    return HIJJET1.kfpj[key2-1][key1-1];
676 }
677
678 //______________________________________________________________________________
679 Float_t THijing::GetPJPX(Int_t key1, Int_t key2) const
680 {
681 // Get values of array PJPX in common HIJJET1
682    if ( key1<1 || key1>300 ) {
683       printf("ERROR in THijing::GetPJPX(key1):\n");
684       printf("      key1=%i is out of range [1..300]\n",key1);
685       return 0;
686    }
687    if ( key2<1 || key2>500 ) {
688       printf("ERROR in THijing::GetPJPX(key1,key2):\n");
689       printf("      key2=%i is out of range [1..500]\n",key2);
690       return 0;
691    }
692    
693    return HIJJET1.pjpx[key2-1][key1-1];
694 }
695
696 //______________________________________________________________________________
697 Float_t THijing::GetPJPY(Int_t key1, Int_t key2) const
698 {
699 // Get values of array PJPY in common HIJJET1
700    if ( key1<1 || key1>300 ) {
701       printf("ERROR in THijing::GetPJPY(key1):\n");
702       printf("      key1=%i is out of range [1..300]\n",key1);
703       return 0;
704    }
705    if ( key2<1 || key2>500 ) {
706       printf("ERROR in THijing::GetPJPY(key1,key2):\n");
707       printf("      key2=%i is out of range [1..500]\n",key2);
708       return 0;
709    }
710    
711    return HIJJET1.pjpy[key2-1][key1-1];
712 }
713
714 //______________________________________________________________________________
715 Float_t THijing::GetPJPZ(Int_t key1, Int_t key2) const
716 {
717 // Get values of array PJPZ in common HIJJET1
718    if ( key1<1 || key1>300 ) {
719       printf("ERROR in THijing::GetPJPZ(key1):\n");
720       printf("      key1=%i is out of range [1..300]\n",key1);
721       return 0;
722    }
723    if ( key2<1 || key2>500 ) {
724       printf("ERROR in THijing::GetPJPZ(key1,key2):\n");
725       printf("      key2=%i is out of range [1..500]\n",key2);
726       return 0;
727    }
728    
729    return HIJJET1.pjpz[key2-1][key1-1];
730 }
731
732 //______________________________________________________________________________
733 Float_t THijing::GetPJPE(Int_t key1, Int_t key2) const
734 {
735 // Get values of array PJPE in common HIJJET1
736    if ( key1<1 || key1>300 ) {
737       printf("ERROR in THijing::GetPJPE(key1):\n");
738       printf("      key1=%i is out of range [1..300]\n",key1);
739       return 0;
740    }
741    if ( key2<1 || key2>500 ) {
742       printf("ERROR in THijing::GetPJPE(key1,key2):\n");
743       printf("      key2=%i is out of range [1..500]\n",key2);
744       return 0;
745    }
746    
747    return HIJJET1.pjpe[key2-1][key1-1];
748 }
749
750 //______________________________________________________________________________
751 Float_t THijing::GetPJPM(Int_t key1, Int_t key2) const
752 {
753 // Get values of array PJPM in common HIJJET1
754    if ( key1<1 || key1>300 ) {
755       printf("ERROR in THijing::GetPJPM(key1):\n");
756       printf("      key1=%i is out of range [1..300]\n",key1);
757       return 0;
758    }
759    if ( key2<1 || key2>500 ) {
760       printf("ERROR in THijing::GetPJPM(key1,key2):\n");
761       printf("      key2=%i is out of range [1..500]\n",key2);
762       return 0;
763    }
764    
765    return HIJJET1.pjpm[key2-1][key1-1];
766 }
767
768 //______________________________________________________________________________
769 Int_t THijing::GetNTJ(Int_t key) const
770 {
771 // Get values of array NTJ in common HIJJET1
772    if ( key<1 || key>300 ) {
773       printf("ERROR in THijing::GetNTJ(key):\n");
774       printf("      key=%i is out of range [1..300]\n",key);
775       return 0;
776    }
777    return HIJJET1.ntj[key-1];
778 }
779
780 //______________________________________________________________________________
781 Int_t THijing::GetKFTJ(Int_t key1, Int_t key2) const
782 {
783 // Get values of array KFTJ in common HIJJET1
784    if ( key1<1 || key1>300 ) {
785       printf("ERROR in THijing::GetKFTJ(key1):\n");
786       printf("      key1=%i is out of range [1..300]\n",key1);
787       return 0;
788    }
789    if ( key2<1 || key2>500 ) {
790       printf("ERROR in THijing::GetKFTJ(key1,key2):\n");
791       printf("      key2=%i is out of range [1..500]\n",key2);
792       return 0;
793    }
794    
795    return HIJJET1.kftj[key2-1][key1-1];
796 }
797
798 //______________________________________________________________________________
799 Float_t THijing::GetPJTX(Int_t key1, Int_t key2) const
800 {
801 // Get values of array PJTX in common HIJJET1
802    if ( key1<1 || key1>300 ) {
803       printf("ERROR in THijing::GetPJTX(key1):\n");
804       printf("      key1=%i is out of range [1..300]\n",key1);
805       return 0;
806    }
807    if ( key2<1 || key2>500 ) {
808       printf("ERROR in THijing::GetPJTX(key1,key2):\n");
809       printf("      key2=%i is out of range [1..500]\n",key2);
810       return 0;
811    }
812    
813    return HIJJET1.pjtx[key2-1][key1-1];
814 }
815
816 //______________________________________________________________________________
817 Float_t THijing::GetPJTY(Int_t key1, Int_t key2) const
818 {
819 // Get values of array PJTY in common HIJJET1
820    if ( key1<1 || key1>300 ) {
821       printf("ERROR in THijing::GetPJTY(key1):\n");
822       printf("      key1=%i is out of range [1..300]\n",key1);
823       return 0;
824    }
825    if ( key2<1 || key2>500 ) {
826       printf("ERROR in THijing::GetPJTY(key1,key2):\n");
827       printf("      key2=%i is out of range [1..500]\n",key2);
828       return 0;
829    }
830    
831    return HIJJET1.pjty[key2-1][key1-1];
832 }
833
834 //______________________________________________________________________________
835 Float_t THijing::GetPJTZ(Int_t key1, Int_t key2) const
836 {
837 // Get values of array PJTZ in common HIJJET1
838    if ( key1<1 || key1>300 ) {
839       printf("ERROR in THijing::GetPJTZ(key1):\n");
840       printf("      key1=%i is out of range [1..300]\n",key1);
841       return 0;
842    }
843    if ( key2<1 || key2>500 ) {
844       printf("ERROR in THijing::GetPJTZ(key1,key2):\n");
845       printf("      key2=%i is out of range [1..500]\n",key2);
846       return 0;
847    }
848    
849    return HIJJET1.pjtz[key2-1][key1-1];
850 }
851
852 //______________________________________________________________________________
853 Float_t THijing::GetPJTE(Int_t key1, Int_t key2) const
854 {
855 // Get values of array PJTE in common HIJJET1
856    if ( key1<1 || key1>300 ) {
857       printf("ERROR in THijing::GetPJTE(key1):\n");
858       printf("      key1=%i is out of range [1..300]\n",key1);
859       return 0;
860    }
861    if ( key2<1 || key2>500 ) {
862       printf("ERROR in THijing::GetPJTE(key1,key2):\n");
863       printf("      key2=%i is out of range [1..500]\n",key2);
864       return 0;
865    }
866    
867    return HIJJET1.pjte[key2-1][key1-1];
868 }
869
870 //______________________________________________________________________________
871 Float_t THijing::GetPJTM(Int_t key1, Int_t key2) const
872 {
873 // Get values of array PJTM in common HIJJET1
874    if ( key1<1 || key1>300 ) {
875       printf("ERROR in THijing::GetPJTM(key1):\n");
876       printf("      key1=%i is out of range [1..300]\n",key1);
877       return 0;
878    }
879    if ( key2<1 || key2>500 ) {
880       printf("ERROR in THijing::GetPJTM(key1,key2):\n");
881       printf("      key2=%i is out of range [1..500]\n",key2);
882       return 0;
883    }
884    
885    return HIJJET1.pjtm[key2-1][key1-1];
886 }
887
888 //====================== access to common HIJJET1 ===============================
889
890 //______________________________________________________________________________
891 Int_t THijing::GetNSG() const
892 {
893 // Get value of NSG in common HIJJET2
894    return HIJJET2.nsg;
895 }
896
897 //______________________________________________________________________________
898 Int_t THijing::GetNJSG(Int_t key) const
899 {
900 // Get values of array NJSG in common HIJJET2
901    if ( key<1 || key>900 ) {
902       printf ("ERROR in THijing:GetNJSG(key): \n ");
903       printf ("      key=%i is out of range [1..900]!\n",key);
904       return 0;
905    }
906
907    return HIJJET2.njsg[key-1];
908
909 }
910
911 //______________________________________________________________________________
912 Int_t THijing::GetIASG(Int_t key1, Int_t key2) const
913 {
914 // Get values of IASG in common HIJJET2
915    if ( key1<1 || key1>900 ) {
916       printf("ERROR in THijing::GetIASG(key1):\n");
917       printf("      key1=%i is out of range [1..900]\n",key1);
918       return 0;
919    }
920    if ( key2<1 || key2>3 ) {
921       printf("ERROR in THijing::GetIASG(key1,key2):\n");
922       printf("      key2=%i is out of range [1..3]\n",key2);
923       return 0;
924    }
925    
926    return HIJJET2.iasg[key2-1][key1-1];
927 }
928
929 //______________________________________________________________________________
930 Int_t THijing::GetK1SG(Int_t key1, Int_t key2) const
931 {
932 // Get values of K1SG in common HIJJET2
933    if ( key1<1 || key1>900 ) {
934       printf("ERROR in THijing::GetK1SG(key1):\n");
935       printf("      key1=%i is out of range [1..900]\n",key1);
936       return 0;
937    }
938    if ( key2<1 || key2>100 ) {
939       printf("ERROR in THijing::GetK1SG(key1,key2):\n");
940       printf("      key2=%i is out of range [1..100]\n",key2);
941       return 0;
942    }
943    
944    return HIJJET2.k1sg[key2-1][key1-1];
945 }
946
947 //______________________________________________________________________________
948 Int_t THijing::GetK2SG(Int_t key1, Int_t key2) const
949 {
950 // Get values of K2SG in common HIJJET2
951    if ( key1<1 || key1>900 ) {
952       printf("ERROR in THijing::GetK2SG(key1):\n");
953       printf("      key1=%i is out of range [1..900]\n",key1);
954       return 0;
955    }
956    if ( key2<1 || key2>100 ) {
957       printf("ERROR in THijing::GetK2SG(key1,key2):\n");
958       printf("      key2=%i is out of range [1..100]\n",key2);
959       return 0;
960    }
961    
962    return HIJJET2.k2sg[key2-1][key1-1];
963 }
964
965 //______________________________________________________________________________
966 Float_t THijing::GetPXSG(Int_t key1, Int_t key2) const
967 {
968 // Get values of PXSG in common HIJJET2
969    if ( key1<1 || key1>900 ) {
970       printf("ERROR in THijing::GetPXSG(key1):\n");
971       printf("      key1=%i is out of range [1..900]\n",key1);
972       return 0;
973    }
974    if ( key2<1 || key2>100 ) {
975       printf("ERROR in THijing::GetPXSG(key1,key2):\n");
976       printf("      key2=%i is out of range [1..100]\n",key2);
977       return 0;
978    }
979    
980    return HIJJET2.pxsg[key2-1][key1-1];
981 }
982
983 //______________________________________________________________________________
984 Float_t THijing::GetPYSG(Int_t key1, Int_t key2) const
985 {
986 // Get values of PYSG in common HIJJET2
987    if ( key1<1 || key1>900 ) {
988       printf("ERROR in THijing::GetPYSG(key1):\n");
989       printf("      key1=%i is out of range [1..900]\n",key1);
990       return 0;
991    }
992    if ( key2<1 || key2>100 ) {
993       printf("ERROR in THijing::GetPYSG(key1,key2):\n");
994       printf("      key2=%i is out of range [1..100]\n",key2);
995       return 0;
996    }
997    
998    return HIJJET2.pysg[key2-1][key1-1];
999 }
1000
1001 //______________________________________________________________________________
1002 Float_t THijing::GetPZSG(Int_t key1, Int_t key2) const
1003 {
1004 // Get values of PZSG in common HIJJET2
1005    if ( key1<1 || key1>900 ) {
1006       printf("ERROR in THijing::GetPZSG(key1):\n");
1007       printf("      key1=%i is out of range [1..900]\n",key1);
1008       return 0;
1009    }
1010    if ( key2<1 || key2>100 ) {
1011       printf("ERROR in THijing::GetPZSG(key1,key2):\n");
1012       printf("      key2=%i is out of range [1..100]\n",key2);
1013       return 0;
1014    }
1015    
1016    return HIJJET2.pzsg[key2-1][key1-1];
1017 }
1018
1019 //______________________________________________________________________________
1020 Float_t THijing::GetPESG(Int_t key1, Int_t key2) const
1021 {
1022 // Get values of PESG in common HIJJET2
1023    if ( key1<1 || key1>900 ) {
1024       printf("ERROR in THijing::GetPESG(key1):\n");
1025       printf("      key1=%i is out of range [1..900]\n",key1);
1026       return 0;
1027    }
1028    if ( key2<1 || key2>100 ) {
1029       printf("ERROR in THijing::GetPESG(key1,key2):\n");
1030       printf("      key2=%i is out of range [1..100]\n",key2);
1031       return 0;
1032    }
1033    
1034    return HIJJET2.pesg[key2-1][key1-1];
1035 }
1036
1037 //______________________________________________________________________________
1038 Float_t THijing::GetPMSG(Int_t key1, Int_t key2) const
1039 {
1040 // Get values of PMSG in common HIJJET2
1041    if ( key1<1 || key1>900 ) {
1042       printf("ERROR in THijing::GetPMSG(key1):\n");
1043       printf("      key1=%i is out of range [1..900]\n",key1);
1044       return 0;
1045    }
1046    if ( key2<1 || key2>100 ) {
1047       printf("ERROR in THijing::GetPMSG(key1,key2):\n");
1048       printf("      key2=%i is out of range [1..100]\n",key2);
1049       return 0;
1050    }
1051    
1052    return HIJJET2.pmsg[key2-1][key1-1];
1053 }
1054
1055 //====================== access to common HISTRNG ===============================
1056
1057 //______________________________________________________________________________
1058 Int_t THijing::GetNFP(Int_t key1, Int_t key2) const
1059 {
1060 // Get values of array NFP in common HISTRNG
1061    if ( key1<1 || key1>300 ) {
1062       printf("ERROR in THijing::GetNFP(key1):\n");
1063       printf("      key1=%i is out of range [1..300]\n",key1);
1064       return 0;
1065    }
1066    if ( key2<1 || key2>15 ) {
1067       printf("ERROR in THijing::GetNFP(key1,key2):\n");
1068       printf("      key2=%i is out of range [1..15]\n",key2);
1069       return 0;
1070    }
1071    
1072    return HISTRNG.nfp[key2-1][key1-1];
1073 }
1074
1075 //______________________________________________________________________________
1076 Float_t THijing::GetPP(Int_t key1, Int_t key2) const
1077 {
1078 // Get values of array PP in common HISTRNG
1079    if ( key1<1 || key1>300 ) {
1080       printf("ERROR in THijing::GetPP(key1):\n");
1081       printf("      key1=%i is out of range [1..300]\n",key1);
1082       return 0;
1083    }
1084    if ( key2<1 || key2>15 ) {
1085       printf("ERROR in THijing::GetPP(key1,key2):\n");
1086       printf("      key2=%i is out of range [1..15]\n",key2);
1087       return 0;
1088    }
1089    
1090    return HISTRNG.pp[key2-1][key1-1];
1091 }
1092
1093 //______________________________________________________________________________
1094 Int_t THijing::GetNFT(Int_t key1, Int_t key2) const
1095 {
1096 // Get values of array NFT in common HISTRNG
1097    if ( key1<1 || key1>300 ) {
1098       printf("ERROR in THijing::GetNFT(key1):\n");
1099       printf("      key1=%i is out of range [1..300]\n",key1);
1100       return 0;
1101    }
1102    if ( key2<1 || key2>15 ) {
1103       printf("ERROR in THijing::GetNFT(key1,key2):\n");
1104       printf("      key2=%i is out of range [1..15]\n",key2);
1105       return 0;
1106    }
1107    
1108    return HISTRNG.nft[key2-1][key1-1];
1109 }
1110
1111 //______________________________________________________________________________
1112 Float_t THijing::GetPT(Int_t key1, Int_t key2) const
1113 {
1114 // Get values of array PT in common HISTRNG
1115    if ( key1<1 || key1>300 ) {
1116       printf("ERROR in THijing::GetPT(key1):\n");
1117       printf("      key1=%i is out of range [1..300]\n",key1);
1118       return 0;
1119    }
1120    if ( key2<1 || key2>15 ) {
1121       printf("ERROR in THijing::GetPT(key1,key2):\n");
1122       printf("      key2=%i is out of range [1..15]\n",key2);
1123       return 0;
1124    }
1125    
1126    return HISTRNG.pt[key2-1][key1-1];
1127 }
1128
1129 void THijing::SetPARJ(Int_t key, Float_t parm) 
1130 {
1131 // Set values of array PARJ in common HISTRNG
1132     if ( key < 1 || key > 200) {
1133         printf("ERROR in THijing::SetPARJ(key,parm):\n");
1134         printf("      key=%i is out of range [1..200]\n",key);
1135     } else {
1136       LUDAT1_HIJING.parj[key-1] = parm;
1137     }
1138 }
1139
1140
1141 void THijing::SetMSTJ(Int_t key, Int_t parm) 
1142 {
1143 // Set values of array MSTJ in common HISTRNG
1144     if ( key < 1 || key > 200) {
1145         printf("ERROR in THijing::SetMSTJ(key,parm):\n");
1146         printf("      key=%i is out of range [1..200]\n",key);
1147     } else {
1148       LUDAT1_HIJING.mstj[key-1] = parm;
1149     }
1150 }
1151
1152 void  THijing::SetMDCY(Int_t key1, Int_t key2, Int_t   parm)
1153 {
1154   // Set value of array MDCY
1155   if ( key1 < 1 || key1 > 500) {
1156     printf("ERROR in THijing::SetMDCY(key1, key2, parm):\n");
1157     printf("      key1=%i is out of range [1..200]\n", key1);
1158   } else if ( key2 < 1 || key2 > 3) {
1159     printf("ERROR in THijing::SetMDCY(key1, key2, parm):\n");
1160     printf("      key2=%i is out of range [1..200]\n", key2);
1161   } else {
1162     LUDAT3_HIJING.mdcy[key2-1][key1-1] = parm;
1163   }
1164   
1165 }
1166
1167 Int_t THijing::GetMDCY(Int_t key1, Int_t key2) 
1168 {
1169   // Get value of array MDCY
1170   if ( key1 < 1 || key1 > 500) {
1171     printf("ERROR in THijing::GetMDCY(key1, key2, parm):\n");
1172     printf("      key1=%i is out of range [1..200]\n", key1);
1173     return -1;
1174   } else if ( key2 < 1 || key2 > 3) {
1175     printf("ERROR in THijing::GetMDCY(key1, key2, parm):\n");
1176     printf("      key2=%i is out of range [1..200]\n", key2);
1177     return -1;
1178   } else {
1179     return (LUDAT3_HIJING.mdcy[key2-1][key1-1]);
1180   }
1181
1182 }
1183 //====================== access to Hijing subroutines =========================
1184
1185
1186 //______________________________________________________________________________
1187 void THijing::Initialize()
1188 {
1189 //////////////////////////////////////////////////////////////////////////////////
1190 // Calls Hijset with the either default parameters or the ones set by the user  //
1191 // via SetEFRM, SetFRAME, SetPROJ, SetTARG, SetIAP, SetIZP, SetIAT, SetIZT      //
1192 //////////////////////////////////////////////////////////////////////////////////
1193
1194    if ( (!strcmp(fFrame.Data(), "CMS     "  )) &&
1195         (!strcmp(fFrame.Data(), "LAB     "  ))){
1196       printf("WARNING! In THijing:Initialize():\n");
1197       printf(" specified frame=%s is neither CMS or LAB\n",fFrame.Data());
1198       printf(" resetting to default \"CMS\" .");
1199       fFrame="CMS";
1200    }
1201
1202    if ( (!strcmp(fProj.Data(), "A       "     )) &&
1203         (!strcmp(fProj.Data(), "P       "     )) &&
1204         (!strcmp(fProj.Data(), "PBAR    "  ))){
1205       printf("WARNING! In THijing:Initialize():\n");
1206       printf(" specified projectile=%s is neither A, P or PBAR\n",fProj.Data());
1207       printf(" resetting to default \"A\" .");
1208       fProj="A";
1209    }
1210
1211    if ( (!strcmp(fTarg.Data(), "A       "     )) &&
1212         (!strcmp(fTarg.Data(), "P       "     )) &&
1213         (!strcmp(fTarg.Data(), "PBAR    "  ))){
1214       printf("WARNING! In THijing:Initialize():\n");
1215       printf(" specified target=%s is neither A, P or PBAR\n",fTarg.Data());
1216       printf(" resetting to default \"A\" .");
1217       fTarg="A";
1218    }
1219
1220    printf(" %s-%s at %f GeV \n",fProj.Data(),fTarg.Data(),fEfrm);
1221
1222    Hijset(fEfrm,fFrame.Data(),fProj.Data(),fTarg.Data(),fIap,fIzp,fIat,fIzt);
1223
1224    printf(" %s-%s at %f GeV \n",fProj.Data(),fTarg.Data(),fEfrm);
1225 }
1226
1227
1228 //______________________________________________________________________________
1229 void THijing::GenerateEvent()
1230 {
1231 // Generates one event;
1232
1233    Hijing(fFrame.Data(),fBmin,fBmax);
1234
1235 }
1236 //______________________________________________________________________________
1237 void THijing::Hijset(float efrm, const char *frame, const char *proj,  
1238                      const char *targ, int iap, int izp, int iat, int izt)
1239 {
1240 // Call HIJING routine HIJSET passing the parameters in a way accepted by
1241 // Fortran routines                                
1242
1243   int s1 = strlen(frame);
1244   int s2 = strlen(proj);
1245   int s3 = strlen(targ);
1246   printf("s1 = %d s2 = %d s3 = %d\n",s1,s2,s3);
1247 #ifndef WIN32 
1248   hijset(efrm, frame, proj, targ, iap, izp, iat, izt, s1, s2, s3);
1249 #else
1250   hijset(efrm, frame, s1, proj, s2, targ, s3, iap, izp, iat, izt);
1251 #endif
1252 }
1253 //______________________________________________________________________________
1254 void THijing::Hijing(const char *frame, float bmin, float bmax)
1255 {
1256 // Call HIJING routine HIJSET passing the parameters in a way accepted by
1257 // Fortran routines                                
1258
1259   int s1 = strlen(frame);
1260   
1261 #ifndef WIN32 
1262   hijing(frame, bmin, bmax, s1);
1263 #else
1264   hijing(frame, s1, bmin, bmax);
1265 #endif
1266 }
1267
1268
1269 Float_t  THijing::Profile(float b)
1270 {
1271 // Call HIJING routine PROFILE 
1272   return profile(b);
1273 }
1274
1275
1276 void  THijing::Rluget(Int_t lfn, Int_t move)
1277 {
1278 // write seed to file
1279   rluget_hijing(lfn, move);
1280 }
1281
1282
1283 void  THijing::Rluset(Int_t lfn, Int_t move)
1284 {
1285 // read seed from file 
1286   rluset_hijing(lfn, move);
1287 }
1288
1289 void  THijing::Pylist(Int_t flag)
1290 {
1291 // call lulist
1292   lulist_hijing(flag);
1293 }