From: morsch Date: Thu, 23 May 2013 15:30:58 +0000 (+0000) Subject: newly implemented model for slow nucleon production X-Git-Url: http://git.uio.no/git/?a=commitdiff_plain;h=58776b75e24c84b77bed993fdfaa649dcaf3b451;p=u%2Fmrichter%2FAliRoot.git newly implemented model for slow nucleon production Chiara Oppedisano --- diff --git a/EVGEN/AliGenSlowNucleons.cxx b/EVGEN/AliGenSlowNucleons.cxx index 7b9a52b7850..7bf44a4b293 100644 --- a/EVGEN/AliGenSlowNucleons.cxx +++ b/EVGEN/AliGenSlowNucleons.cxx @@ -140,7 +140,7 @@ void AliGenSlowNucleons::FinishRun() TCanvas *c = new TCanvas("c","Canvas 1",400,10,600,700); c->Divide(2,1); c->cd(1); - fDebugHist1->Draw(); + fDebugHist1->Draw("colz"); c->cd(2); fDebugHist2->Draw(); c->cd(3); @@ -164,11 +164,12 @@ void AliGenSlowNucleons::Generate() // Int_t nnw = fCollisionGeometry->NNw(); // Int_t nwnw = fCollisionGeometry->NwNw(); - fSlowNucleonModel->GetNumberOfSlowNucleons(fCollisionGeometry, fNgp, fNgn, fNbp, fNbn); + //fSlowNucleonModel->GetNumberOfSlowNucleons(fCollisionGeometry, fNgp, fNgn, fNbp, fNbn); + fSlowNucleonModel->GetNumberOfSlowNucleons2(fCollisionGeometry, fNgp, fNgn, fNbp, fNbn); if (fDebug) { //printf("Collision Geometry %f %d %d %d %d\n", b, nn, nwn, nnw, nwnw); printf("Slow nucleons: %d grayp %d grayn %d blackp %d blackn \n", fNgp, fNgn, fNbp, fNbn); - fDebugHist1->Fill(Float_t(fNgp + fNgn + fNbp + fNbn), fCollisionGeometry->NwN(), 1.); + fDebugHist1->Fill(Float_t(fNgp + fNgn + fNbp + fNbn), fCollisionGeometry->NNw(), 1.); fDebugHist2->Fill(Float_t(fNgp + fNgn + fNbp + fNbn), b, 1.); } } @@ -195,7 +196,7 @@ void AliGenSlowNucleons::Generate() GenerateSlow(fCharge, fTemperatureG, fBetaSourceG, p, theta); if (fDebug) fCosThetaGrayHist->Fill(TMath::Cos(theta)); PushTrack(fTrackIt, -1, kf, p, origin, polar, - time, kPNoProcess, nt, 1.,1); + time, kPNoProcess, nt, 1.,-2); KeepTrack(nt); } // @@ -207,7 +208,7 @@ void AliGenSlowNucleons::Generate() GenerateSlow(fCharge, fTemperatureG, fBetaSourceG, p, theta); if (fDebug) fCosThetaGrayHist->Fill(TMath::Cos(theta)); PushTrack(fTrackIt, -1, kf, p, origin, polar, - time, kPNoProcess, nt, 1.,1); + time, kPNoProcess, nt, 1.,-2); KeepTrack(nt); } // @@ -218,7 +219,7 @@ void AliGenSlowNucleons::Generate() for(i = 0; i < fNbp; i++) { GenerateSlow(fCharge, fTemperatureB, fBetaSourceB, p, theta); PushTrack(fTrackIt, -1, kf, p, origin, polar, - time, kPNoProcess, nt, 1.,1); + time, kPNoProcess, nt, 1.,-1); KeepTrack(nt); } // @@ -229,7 +230,7 @@ void AliGenSlowNucleons::Generate() for(i = 0; i < fNbn; i++) { GenerateSlow(fCharge, fTemperatureB, fBetaSourceB, p, theta); PushTrack(fTrackIt, -1, kf, p, origin, polar, - time, kPNoProcess, nt, 1.,1); + time, kPNoProcess, nt, 1.,-1); KeepTrack(nt); } } diff --git a/EVGEN/AliSlowNucleonModel.h b/EVGEN/AliSlowNucleonModel.h index 38b85a3cfc5..d78d7531ec9 100644 --- a/EVGEN/AliSlowNucleonModel.h +++ b/EVGEN/AliSlowNucleonModel.h @@ -16,6 +16,9 @@ public: virtual void GetNumberOfSlowNucleons(AliCollisionGeometry* /*geo*/, Int_t& /*ngp*/, Int_t& /*ngn*/, Int_t& /*nbp*/, Int_t& /*nbn*/) const {;} + virtual void GetNumberOfSlowNucleons2(AliCollisionGeometry* /*geo*/, + Int_t& /*ngp*/, Int_t& /*ngn*/, + Int_t& /*nbp*/, Int_t& /*nbn*/) const {;} protected: ClassDef(AliSlowNucleonModel,1) // Gray Particle Model diff --git a/EVGEN/AliSlowNucleonModelExp.cxx b/EVGEN/AliSlowNucleonModelExp.cxx index 9598cb2bf97..5b73473d122 100644 --- a/EVGEN/AliSlowNucleonModelExp.cxx +++ b/EVGEN/AliSlowNucleonModelExp.cxx @@ -26,13 +26,14 @@ #include "AliSlowNucleonModelExp.h" #include "AliCollisionGeometry.h" #include +#include ClassImp(AliSlowNucleonModelExp) AliSlowNucleonModelExp::AliSlowNucleonModelExp(): fP(82), - fN (208 - 82), + fN (126), fAlphaGray(2.3), fAlphaBlack(3.6), fApplySaturation(kTRUE), @@ -45,7 +46,7 @@ AliSlowNucleonModelExp::AliSlowNucleonModelExp(): // printf("\n\nInitializing slow nucleon model with parameters:\n"); printf(" \t alpha_{gray} %1.2f alpha_{black} %1.2f\n",fAlphaGray, fAlphaBlack); - printf(" \t SATURATION %d w. %d (gray) %d (black) \n\n",fApplySaturation,fnGraySaturation,fnBlackSaturation); + //printf(" \t SATURATION %d w. %d (gray) %d (black) \n\n",fApplySaturation,fnGraySaturation,fnBlackSaturation); } @@ -57,7 +58,7 @@ void AliSlowNucleonModelExp::GetNumberOfSlowNucleons(AliCollisionGeometry* geo, // // Number of collisions - Int_t nu = geo->NN() + geo->NwN() + geo->NNw(); + Float_t nu = geo->NN() + geo->NwN() + geo->NNw(); // Mean number of gray nucleons @@ -90,7 +91,79 @@ void AliSlowNucleonModelExp::GetNumberOfSlowNucleons(AliCollisionGeometry* geo, // black protons p = nBlackProtons/fP; nbp = gRandom->Binomial((Int_t) fP, p); - + +} + +void AliSlowNucleonModelExp::GetNumberOfSlowNucleons2(AliCollisionGeometry* geo, + Int_t& ngp, Int_t& ngn, Int_t & nbp, Int_t & nbn) const +{ +// +// Return the number of black and gray nucleons +// +// Number of collisions + + // based on E910 model ================================================================ + + Float_t nu = (Float_t) (geo->NN() + geo->NwN() + geo->NNw()); + // + nu = nu+1.*gRandom->Rndm(); + // + Float_t poverpd = 0.843; + Float_t zAu2zPb = 82./79.; + Float_t nGrayp = (-0.27 + 0.63 * nu - 0.0008 *nu *nu)*poverpd*zAu2zPb; + +// gray protons + Double_t p; + p = nGrayp/fP; + ngp = gRandom->Binomial((Int_t) fP, p); + //ngp = gRandom->Gaus(nGrayp, TMath::Sqrt(fP*p*(1-p))); + if(nGrayp<0.) ngp=0; + + //Float_t blackovergray = 3./7.;// from spallation + Float_t blackovergray = 0.65; // from COSY + Float_t nBlackp = blackovergray*nGrayp; + +// black protons + p = nBlackp/fP; + nbp = gRandom->Binomial((Int_t) fP, p); + //nbp = gRandom->Gaus(nBlackp, TMath::Sqrt(fP*p*(1-p))); + if(nBlackp<0.) nbp=0; + + if(nu<3.){ + nGrayp = -0.836 + 0.9112 *nu - 0.05381 *nu *nu; + nBlackp = blackovergray*nGrayp; + } + + Float_t nGrayNeutrons = 0.; + Float_t nBlackNeutrons = 0.; + Float_t cp = (nGrayp+nBlackp)/0.24; + + if(cp>0.){ + Float_t nSlow = 51.5+469.2/(-8.762-cp); + //if(cp<2.5) nSlow = 1+(9.9-1)/(2.5-0)*(cp-0); + if(cp<3.) nSlow = 0.+(11.6-0.)/(3.-0.)*(cp-0.); + + nGrayNeutrons = nSlow * 0.1; + nBlackNeutrons = nSlow - nGrayNeutrons; + } + else{ + // Sikler "pasturato" + nGrayNeutrons = 0.47 * fAlphaGray * nu; + nBlackNeutrons = 0.88 * fAlphaBlack * nu; + printf("nslowp=0 -> ncoll = %1.0f -> ngrayn = %1.0f nblackn = %1.0f \n", nu, nGrayNeutrons, nBlackNeutrons); + } + +// gray neutrons + p = nGrayNeutrons/fN; +// ngn = gRandom->Binomial((Int_t) fN, p); + ngn = gRandom->Gaus(nGrayNeutrons, TMath::Sqrt(fN*p*(1-p))); + +// black neutrons + p = nBlackNeutrons/fN; +// nbn = gRandom->Binomial((Int_t) fN, p); + nbn = gRandom->Gaus(nBlackNeutrons, TMath::Sqrt(fN*p*(1-p))); + + } void AliSlowNucleonModelExp::SetParameters(Float_t alpha1, Float_t alpha2) diff --git a/EVGEN/AliSlowNucleonModelExp.h b/EVGEN/AliSlowNucleonModelExp.h index 83a240dfa64..4755e44c4a7 100644 --- a/EVGEN/AliSlowNucleonModelExp.h +++ b/EVGEN/AliSlowNucleonModelExp.h @@ -22,6 +22,8 @@ class AliSlowNucleonModelExp : public AliSlowNucleonModel virtual ~AliSlowNucleonModelExp(){;} virtual void GetNumberOfSlowNucleons(AliCollisionGeometry* geo, Int_t& ngp, Int_t& ngn, Int_t& nbp, Int_t& nbn) const; + virtual void GetNumberOfSlowNucleons2(AliCollisionGeometry* geo, + Int_t& ngp, Int_t& ngn, Int_t& nbp, Int_t& nbn) const; virtual void SetParameters(Float_t alpha1, Float_t alpha2); virtual void SetSaturation(Bool_t saturation) {fApplySaturation = saturation;} virtual void SetSaturationParams(Int_t ngray=15, Int_t nblack=28)