THijing/THijing.cxx

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