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