Chiara O. implemented a way to bypass HIJING internal calculation
authormorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Sat, 26 Jan 2013 08:42:30 +0000 (08:42 +0000)
committermorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Sat, 26 Jan 2013 08:42:30 +0000 (08:42 +0000)
of the sigmaNN so that one can set the value from outside, and to switch off the elastic collisions.
By default everything should be as before.

Constantin L.  then improved the calculation of the geometrical Npart,
which now is very close to the MC Glauber one. It still differs
because HIJING normalizes internally to collisions with the eikonal
approach for which the geometrical Npart can be zero (if we would
accept events for which the geometrical Npart > 0 then the result would be equal to the one of the MC Glauber).

HIJING/hijing1_36/hijing.F
HIJING/hijing1_36/hijset.F
THijing/AliGenHijing.cxx
THijing/AliGenHijing.h

index 04659cf..0040018 100644 (file)
@@ -154,7 +154,7 @@ C****************************************************************
        SUBROUTINE HIJING(FRAME,BMIN0,BMAX0)
        CHARACTER FRAME*8
        DIMENSION SCIP(300,300),RNIP(300,300),SJIP(300,300),JTP(3),
-     &                 IPCOL(90000),ITCOL(90000)
+     &                 IPCOL(90000),ITCOL(90000),SCIP2(300,300)
 #define BLANKET_SAVE
 #include "hiparnt.inc"
 C
@@ -342,6 +342,7 @@ C
           SCIP(JP,JT)=-1.0
           B2=(YP(1,JP)+BBX-YT(1,JT))**2+(YP(2,JP)+BBY-YT(2,JT))**2
           R2=B2*HIPR1(40)/HIPR1(31)/0.1
+          SCIP2(JP,JT)=R2
 C              ********mb=0.1*fm, YP is in fm,HIPR1(31) is in mb
           RRB1=MIN((YP(1,JP)**2+YP(2,JP)**2)
      &          /1.2**2/REAL(IHNT2(1))**0.6666667,1.0)
@@ -376,6 +377,9 @@ C                   ********perform elastic collisions
           RNIP(JP,JT)=RANTOT
           SJIP(JP,JT)=HINT1(18)
           NCOLT=NCOLT+1
+           if (R2.GT.2.D0) THEN
+              write (8,*) R2
+           ENDIF
           IPCOL(NCOLT)=JP
           ITCOL(NCOLT)=JT
 70     CONTINUE
@@ -399,7 +403,7 @@ c *** cl glauber ***
               xmeana=xmeana+YP(1,JP)
               ymeana=ymeana+YP(2,JP)
               DO 1120 JT=1,IHNT2(3)
-                 IF(SCIP(JP,JT).GT.-1.0D0) THEN
+                 IF(SCIP2(JP,JT).LT.2.0D0) THEN
                     npart=npart+1
                     xmeanp=xmeanp+YP(1,JP)
                     ymeanp=ymeanp+YP(2,JP)
@@ -415,7 +419,7 @@ c *** cl glauber ***
               xmeanb=xmeanb+YT(1,JT)
               ymeanb=ymeanb+YT(2,JT)
               DO 1140 JP=1,IHNT2(1)
-                 IF(SCIP(JP,JT).GT.-1.0D0) THEN
+                 IF(SCIP2(JP,JT).LT.2.0D0) THEN
                     npart=npart+1
                     xmeanp=xmeanp+YT(1,JT)
                     ymeanp=ymeanp+YT(2,JT)
@@ -457,6 +461,7 @@ c *** cl glauber ***
 
 C              ********total number interactions proj and targ has
 C                              suffered
+
        IF(NCOLT.EQ.0) THEN
           NLOP=NLOP+1
            IF(NLOP.LE.20.OR.
index 38c8651..65f861d 100644 (file)
@@ -188,8 +188,8 @@ C
           HIDAT(J)=HIDAT0(J,I-1)+(HIDAT0(J,I)-HIDAT0(J,I-1))
      &    *(HINT1(1)-HIDAT0(10,I-1))/(HIDAT0(10,I)-HIDAT0(10,I-1))
 40     CONTINUE
-       HIPR1(31)=HIDAT(5)
-       HIPR1(30)=2.0*HIDAT(5)
+       IF(HIPR1(31).EQ.0) HIPR1(31)=HIDAT(5)
+       HIPR1(30)=2.0*HIPR1(31)
 C
 C
        CALL HIJCRS
@@ -218,6 +218,7 @@ C           ********booking for x distribution of valence quarks
 
        write (*,*),'Average jet cross section (mb): " ', HINT1(10)
        write (*,*),'Jet cross section (mb): " ', HINT1(11)
+       write (*,*),'Geometrical N size (mb): ', HIPR1(31)
        write (*,*),'Inelastic NN cross section (mb): ', HINT1(12)
         write (*,*),'Total NN cross section (mb): ',HINT1(13) 
        RETURN
index df3ccea..43cbbd8 100644 (file)
@@ -74,7 +74,9 @@ AliGenHijing::AliGenHijing()
      fLHC(kFALSE),
      fRandomPz(kFALSE),
      fNoHeavyQuarks(kFALSE),
-     fHeader(AliGenHijingEventHeader("Hijing"))
+     fHeader(AliGenHijingEventHeader("Hijing")),
+     fSigmaNN(-1),
+     fNoElas(0)
 {
   // Constructor
   fEnergyCMS = 5500.;
@@ -118,7 +120,9 @@ AliGenHijing::AliGenHijing(Int_t npart)
      fLHC(kFALSE),
      fRandomPz(kFALSE),
      fNoHeavyQuarks(kFALSE),
-     fHeader(AliGenHijingEventHeader("Hijing"))
+     fHeader(AliGenHijingEventHeader("Hijing")),
+     fSigmaNN(-1),
+     fNoElas(0)
 {
 // Default PbPb collisions at 5. 5 TeV
 //
@@ -158,8 +162,14 @@ void AliGenHijing::Init()
     fHijing->SetIHPR2(21, fKeep);
     fHijing->SetHIPR1(8,  fPtHardMin);         
     fHijing->SetHIPR1(9,  fPtHardMax);         
-    fHijing->SetHIPR1(10, fPtMinJet);  
+    fHijing->SetHIPR1(10, fPtMinJet); 
+    if (fSigmaNN>0)
+      fHijing->SetHIPR1(31, fSigmaNN/2.);      
     fHijing->SetHIPR1(50, fSimpleJet);
+    //
+    // Switching off elastic scattering
+    if (fNoElas)
+      fHijing->SetIHPR2(14, 0);
 //
 //  Quenching
 //
index 010ae41..0dd4366 100644 (file)
@@ -56,11 +56,12 @@ class AliGenHijing : public AliGenMC
        {fEtaMinJet = etamin; fEtaMaxJet = etamax;}
     virtual void    SetJetPhiRange(Float_t phimin = -180., Float_t phimax = 180.)
        {fPhiMinJet = TMath::Pi()*phimin/180.; fPhiMaxJet = TMath::Pi()*phimax/180.;}
-    virtual void    SetBoostLHC(Int_t flag = 0)         {fLHC        = flag;}
-    virtual void    SetRandomPz(Bool_t flag = 0)        {fRandomPz   = flag;}
+    virtual void    SetBoostLHC(Int_t flag = 0)       {fLHC        = flag;}
+    virtual void    SetRandomPz(Bool_t flag = 0)      {fRandomPz   = flag;}
     virtual void    SwitchOffHeavyQuarks(Bool_t flag = kTRUE) {fNoHeavyQuarks = flag;}
-    
-           
+    virtual void    SetSigmaNN(Float_t val)           {fSigmaNN    = val;}    
+    virtual void    SetNoElas(Bool_t b)               {fNoElas     = b; }          
+
 // Getters
     virtual TString GetReferenceFrame()  const {return fFrame;}
     virtual void    GetImpactParameterRange(Float_t& bmin, Float_t& bmax) const
@@ -79,6 +80,7 @@ class AliGenHijing : public AliGenMC
     virtual void    GetJetPhiRange(Float_t& phimin, Float_t& phimax)      const
        {phimin = fPhiMinJet*180./TMath::Pi(); phimax = fPhiMaxJet*180./TMath::Pi();}
      THijing       *GetTHijing()                         const {return fHijing;}
+    virtual Float_t GetSigmaNN()                         const {return fSigmaNN;}
 
 // Physics Routines
     virtual Bool_t  ProvidesCollisionGeometry() const {return kTRUE;}
@@ -128,7 +130,9 @@ class AliGenHijing : public AliGenMC
     Bool_t      fRandomPz;       // Randomise sign of pz  event by event
     Bool_t      fNoHeavyQuarks;  // If true no heavy quarks are produced
     AliGenHijingEventHeader     fHeader; // MC Header
-    
+    Float_t     fSigmaNN;        // If not -1 set sigmaNN (HIPR1) 
+    Bool_t      fNoElas;         // If true switch off elastic scattering
+
  private:
     AliGenHijing(const AliGenHijing &Hijing);
     AliGenHijing &  operator=(const AliGenHijing & rhs);
@@ -140,11 +144,6 @@ class AliGenHijing : public AliGenMC
     // check if stable
     Bool_t Stable(const TParticle*  particle) const;
     
-    ClassDef(AliGenHijing, 8) // AliGenerator interface to Hijing
+    ClassDef(AliGenHijing, 9) // AliGenerator interface to Hijing
 };
 #endif
-
-
-
-
-