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