EVGEN/AliSlowNucleonModelExp.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /* $Id$ */
  17
  18 //
  19 // Experimental data inspired Gray Particle Model for p-Pb collisions
  20 // The number of gray nucleons  is proportional to the number of collisions.
  21 // The number of black nucleons is proportional to the number of collisions
  22 // Fluctuations are calculated from a binomial distribution.
  23 // Author: A.Morsch
  24 //
  25
  26 #include "AliSlowNucleonModelExp.h"
  27 #include "AliCollisionGeometry.h"
  28 #include <TRandom.h>
  29 #include <TMath.h>
  30
  31 ClassImp(AliSlowNucleonModelExp)
  32
  33
  34 AliSlowNucleonModelExp::AliSlowNucleonModelExp():
  35     fP(82),
  36     fN (126),
  37     fAlphaGray(2.3),
  38     fAlphaBlack(3.6),
  39     fApplySaturation(kTRUE),
  40     fnGraySaturation(15),
  41     fnBlackSaturation(28),
  42     fLCPparam(0.585),
  43     fSigmaSmear(0.25)
  44 {
  45   //
  46   // Default constructor
  47   //
  48   //
  49   fSlownparam[0] = 60.;
  50   fSlownparam[1] = 469.2;
  51   fSlownparam[2] = 8.762;
  52   /*printf("\n\n ******** Initializing slow nucleon model with parameters:\n");
  53   printf(" \t alpha_{gray} %1.2f  alpha_{black} %1.2f\n",fAlphaGray, fAlphaBlack);
  54   printf(" \t SATURATION %d w. %d (gray) %d (black) \n\n",fApplySaturation,fnGraySaturation,fnBlackSaturation);
  55   printf(" \t LCP parameter %f   Slown parameters = {%f, %f,
  56   %f}\n\n",fLCPparam,fSlownparam[0],fSlownparam[1],fSlownparam[2]); */
  57 }
  58
  59
  60 void AliSlowNucleonModelExp::GetNumberOfSlowNucleons(AliCollisionGeometry* geo,
  61                                                       Int_t& ngp, Int_t& ngn, Int_t & nbp, Int_t & nbn) const
  62 {
  63 //
  64 // Return the number of black and gray nucleons
  65 //
  66 // Number of collisions
  67
  68     Float_t nu = geo->NN() + geo->NwN() + geo->NNw();
  69
  70 // Mean number of gray nucleons
  71
  72     Float_t nGray         = fAlphaGray * nu;
  73     Float_t nGrayNeutrons = nGray * fN / (fN + fP);
  74     Float_t nGrayProtons  = nGray - nGrayNeutrons;
  75
  76 // Mean number of black nucleons
  77     Float_t nBlack  = 0.;
  78     if(!fApplySaturation || (fApplySaturation && nGray<fnGraySaturation)) nBlack = fAlphaBlack * nu;
  79     else if(fApplySaturation && nGray>=fnGraySaturation) nBlack = fnBlackSaturation;
  80     Float_t nBlackNeutrons = nBlack * 0.84;
  81     Float_t nBlackProtons  = nBlack - nBlackNeutrons;
  82
  83 // Actual number (including fluctuations) from binomial distribution
  84     Double_t p;
  85
  86 //  gray neutrons
  87     p =  nGrayNeutrons/fN;
  88     ngn = gRandom->Binomial((Int_t) fN, p);
  89
  90 //  gray protons
  91     p =  nGrayProtons/fP;
  92     ngp = gRandom->Binomial((Int_t) fP, p);
  93
  94 //  black neutrons
  95     p =  nBlackNeutrons/fN;
  96     nbn = gRandom->Binomial((Int_t) fN, p);
  97
  98 //  black protons
  99     p =  nBlackProtons/fP;
 100     nbp = gRandom->Binomial((Int_t) fP, p);
 101
 102 }
 103
 104 void AliSlowNucleonModelExp::GetNumberOfSlowNucleons2(AliCollisionGeometry* geo,
 105                                                       Int_t& ngp, Int_t& ngn, Int_t & nbp, Int_t & nbn) const
 106 {
 107 //
 108 // Return the number of black and gray nucleons
 109 //
 110 // Number of collisions
 111
 112    // based on E910 model ================================================================
 113
 114    Float_t nu = (Float_t) (geo->NN() + geo->NwN() + geo->NNw());
 115    //
 116    //nu = nu+1.*gRandom->Rndm();
 117    nu = gRandom->Gaus(nu, 0.5);
 118    if(nu<0.) nu=0.;
 119    //
 120    Float_t  poverpd = 0.843;
 121    Float_t  zAu2zPb = 82./79.;
 122    Float_t  nGrayp = (-0.27 + 0.63 * nu - 0.0008 *nu *nu)*poverpd*zAu2zPb;
 123
 124 //  gray protons
 125     Double_t p;
 126     p =  nGrayp/fP;
 127     ngp = gRandom->Binomial((Int_t) fP, p);
 128     //ngp = gRandom->Gaus(nGrayp, TMath::Sqrt(fP*p*(1-p)));
 129     if(nGrayp<0.) ngp=0;
 130
 131     //Float_t blackovergray = 3./7.;// from spallation
 132     Float_t blackovergray = 0.65; // from COSY
 133     Float_t nBlackp  = blackovergray*nGrayp;
 134
 135 //  black protons
 136     p =  nBlackp/fP;
 137     nbp = gRandom->Binomial((Int_t) fP, p);
 138     //nbp = gRandom->Gaus(nBlackp, TMath::Sqrt(fP*p*(1-p)));
 139     if(nBlackp<0.) nbp=0;
 140
 141     if(nu<3.){
 142       nGrayp = -0.836 + 0.9112 *nu - 0.05381 *nu *nu;
 143       nBlackp  = blackovergray*nGrayp;
 144     }
 145
 146     //printf(" \t Using LCP parameter %f   Slown parameters = {%f, %f, %f}\n\n",fLCPparam,fSlownparam[0],fSlownparam[1],fSlownparam[2]);
 147     Float_t nGrayNeutrons = 0.;
 148     Float_t nBlackNeutrons = 0.;
 149     Float_t cp = (nGrayp+nBlackp)/fLCPparam;
 150
 151     if(cp>0.){
 152       Float_t nSlow      = fSlownparam[0]+fSlownparam[1]/(-fSlownparam[2]-cp);
 153       Float_t paramRetta = fSlownparam[0]+fSlownparam[1]/(-fSlownparam[2]-3);
 154       if(cp<3.) nSlow = 0.+(paramRetta-0.)/(3.-0.)*(cp-0.);
 155
 156       nGrayNeutrons = nSlow * 0.1;
 157       nBlackNeutrons = nSlow - nGrayNeutrons;
 158     }
 159     else{
 160       // Sikler "pasturato" (qui non entra mai!!!!)
 161       nGrayNeutrons = 0.47 * fAlphaGray *  nu;
 162       nBlackNeutrons = 0.88 * fAlphaBlack * nu;
 163       //printf("nslowp=0 -> ncoll = %1.0f -> ngrayn = %1.0f  nblackn = %1.0f \n", nu, nGrayNeutrons, nBlackNeutrons);
 164     }
 165
 166 //  gray neutrons
 167     p =  nGrayNeutrons/fN;
 168 //    ngn = gRandom->Binomial((Int_t) fN, p);
 169     ngn = gRandom->Gaus(nGrayNeutrons, TMath::Sqrt(fN*p*(1-p)));
 170
 171 //  black neutrons
 172     p =  nBlackNeutrons/fN;
 173 //    nbn = gRandom->Binomial((Int_t) fN, p);
 174     nbn = gRandom->Gaus(nBlackNeutrons, TMath::Sqrt(fN*p*(1-p)));
 175
 176
 177 }
 178
 179 void AliSlowNucleonModelExp::GetNumberOfSlowNucleons2s(AliCollisionGeometry* geo,
 180                                                       Int_t& ngp, Int_t& ngn, Int_t & nbp, Int_t & nbn) const
 181 {
 182 //
 183 // Return the number of black and gray nucleons
 184 //
 185 // Number of collisions
 186
 187    // based on E910 model ================================================================
 188
 189    Float_t nu = (Float_t) (geo->NN() + geo->NwN() + geo->NNw());
 190    //
 191    Float_t  poverpd = 0.843;
 192    Float_t  zAu2zPb = 82./79.;
 193    Float_t  grayp = (-0.27 + 0.63 * nu - 0.0008 *nu *nu)*poverpd*zAu2zPb;
 194    Float_t  nGrayp = gRandom->Gaus(grayp, fSigmaSmear);
 195    if(nGrayp<0.) nGrayp=0.;
 196
 197 //  gray protons
 198     Double_t p=0.;
 199     p = nGrayp/fP;
 200     ngp = gRandom->Binomial((Int_t) fP, p);
 201     //ngp = gRandom->Gaus(nGrayp, TMath::Sqrt(fP*p*(1-p)));
 202     if(nGrayp<0.) ngp=0;
 203
 204     //Float_t blackovergray = 3./7.;// from spallation
 205     Float_t blackovergray = 0.65; // from COSY
 206     //Float_t blackp  = blackovergray*grayp;
 207     //Float_t nBlackp = gRandom->Gaus(nblackp, fSigmaSmear);
 208     Float_t nBlackp = blackovergray*nGrayp;
 209     if(nBlackp<0.) nBlackp=0.;
 210
 211 //  black protons
 212     p =  nBlackp/fP;
 213     nbp = gRandom->Binomial((Int_t) fP, p);
 214     //nbp = gRandom->Gaus(nBlackp, TMath::Sqrt(fP*p*(1-p)));
 215     if(nBlackp<0.) nbp=0;
 216
 217     Float_t nGrayNeutrons = 0.;
 218     Float_t nBlackNeutrons = 0.;
 219     Float_t cp = (nGrayp+nBlackp)/fLCPparam;
 220
 221     if(cp>0.){
 222       Float_t nSlow = fSlownparam[0]+fSlownparam[1]/(-fSlownparam[2]-cp);
 223
 224       nGrayNeutrons = nSlow * 0.1;
 225       nBlackNeutrons = nSlow - nGrayNeutrons;
 226     }
 227     else{
 228       // Sikler "pasturato" (qui non entra mai!!!!)
 229       nGrayNeutrons = 0.47 * fAlphaGray *  nu;
 230       nBlackNeutrons = 0.88 * fAlphaBlack * nu;
 231       //printf("nslowp=0 -> ncoll = %1.0f -> ngrayn = %1.0f  nblackn = %1.0f \n", nu, nGrayNeutrons, nBlackNeutrons);
 232     }
 233     //
 234     if(nGrayNeutrons<0.) nGrayNeutrons=0.;
 235     if(nBlackNeutrons<0.) nBlackNeutrons=0.;
 236
 237 //  gray neutrons
 238     p =  nGrayNeutrons/fN;
 239 //    ngn = gRandom->Binomial((Int_t) fN, p);
 240     ngn = gRandom->Gaus(nGrayNeutrons, TMath::Sqrt(fN*p*(1-p)));
 241     if(nGrayNeutrons<0.) ngn=0;
 242
 243 //  black neutrons
 244     p =  nBlackNeutrons/fN;
 245 //    nbn = gRandom->Binomial((Int_t) fN, p);
 246     nbn = gRandom->Gaus(nBlackNeutrons, TMath::Sqrt(fN*p*(1-p)));
 247     if(nBlackNeutrons<0.) nbn=0;
 248
 249 }
 250
 251 void AliSlowNucleonModelExp::SetParameters(Float_t alpha1, Float_t alpha2)
 252 {
 253     // Set the model parameters
 254     fAlphaGray  = alpha1;
 255     fAlphaBlack = alpha2;
 256 }
 257