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