]>
Commit | Line | Data |
---|---|---|
e0b1fdb4 | 1 | void Fragment(Int_t nev=0, Int_t debug=0) |
2 | { | |
3 | gROOT->Reset(); | |
e0b1fdb4 | 4 | |
5 | Float_t b; | |
6 | Int_t nspectp, nspectn, nspectpfree, nspectnfree; | |
7 | Int_t zz[100], nn[100], nalpha, ztot, ntot; | |
8 | void spectator(Float_t, Int_t*, Int_t*); | |
9 | ||
10 | TH2F *hsp = new TH2F("hsp","Spectator protons vs b",100,0.,20.,100,0.,150.); | |
11 | hsp -> SetXTitle("b (fm)"); | |
12 | hsp -> SetYTitle("Num. of spectator protons"); | |
13 | ||
14 | TH2F *hsn = new TH2F("hsn","Spectator neutrons vs b",100,0.,20.,100,0.,150.); | |
15 | hsn -> SetXTitle("b (fm)"); | |
16 | hsn -> SetYTitle("Num. of spectator neutrons"); | |
17 | ||
18 | TH2F *hFragp = new TH2F("hFragp","Free spectator protons vs b",100,0.,20.,100,0.,150.); | |
19 | hFragp -> SetXTitle("b (fm)"); | |
20 | hFragp -> SetYTitle("Num. of free spectator protons"); | |
21 | ||
22 | TH2F *hFragn = new TH2F("hFragn","Free spectator neutrons vs b",100,0.,20.,100,0.,150.); | |
23 | hFragn -> SetXTitle("b (fm)"); | |
24 | hFragn -> SetYTitle("Num. of free spectator neutrons"); | |
25 | ||
26 | TF1 *funb = new TF1("funb","x",0,15); | |
27 | for(Int_t ievent=0; ievent<nev; ievent++){ | |
28 | if(ievent%1000==0){printf("Processing event %d\n",ievent);} | |
29 | b= Float_t(funb->GetRandom()); | |
30 | spectator(b, &nspectp, &nspectn); | |
31 | hsp -> Fill(b,nspectp); | |
32 | hsn -> Fill(b,nspectn); | |
33 | AliZDCFragment *gallio = new AliZDCFragment(b); | |
34 | for(Int_t j=0; j<=99; j++){ | |
35 | zz[j] =0; | |
36 | nn[j] =0; | |
37 | } | |
38 | // | |
39 | // Generation of fragments | |
25bd4833 | 40 | gallio->GenerateIMF(); |
e0b1fdb4 | 41 | Int_t NFrag = gallio->GetFragmentNum(); |
42 | // | |
43 | // Attach neutrons | |
44 | ztot=0; | |
45 | ntot=0; | |
25bd4833 | 46 | gallio->AttachNeutrons(); |
e0b1fdb4 | 47 | nspectpfree = nspectp-ztot-2*nalpha; |
48 | nspectnfree = nspectn-ntot-2*nalpha; | |
49 | hFragp -> Fill(b,nspectpfree); | |
50 | hFragn -> Fill(b,nspectnfree); | |
51 | // | |
52 | ||
53 | if(debug == 1){ | |
54 | printf("\n b = %f\n",b); | |
55 | printf("\n #spect p = %d, #spect n = %d\n",nspectp,nspectn); | |
56 | printf("\n #spect p free = %d, #spect n free = %d\n",nspectpfree,nspectnfree); | |
57 | printf("\n #fragments = %f\n",NFrag); | |
58 | for(Int_t i=0; i<NFrag; i++){ | |
59 | printf("\n ZZ[%d] = %d, NN[%d] = %d\n",i,zz[i],i,nn[i]); | |
60 | } | |
61 | printf("\n NAlpha = %d, Ztot = %d, Ntot = %d\n\n",nalpha,ztot,ntot); | |
62 | } | |
63 | delete gallio; | |
64 | } //Event loop | |
65 | ||
66 | TCanvas *c1 = new TCanvas("c1","Fragmentation",0,10,580,700); | |
67 | c1->cd(); | |
68 | c1->SetFillColor(29); | |
69 | c1->Divide(2,2); | |
70 | c1->cd(1); | |
71 | hsp -> Draw("box"); | |
72 | c1->cd(2); | |
73 | hsn -> Draw("box"); | |
74 | c1->cd(3); | |
75 | hFragp -> Draw("box"); | |
76 | c1->cd(4); | |
77 | hFragn -> Draw("box"); | |
78 | ||
79 | } | |
80 | ||
81 | ||
82 | void spectator(Float_t b, Int_t* NSpectp, Int_t* NSpectn) | |
83 | { | |
84 | Float_t SppvsB[6] = {3.633,-1.518,1.360,-.4899e-1,-.2398e-2,.1066e-3}; | |
85 | Float_t SpnvsB[6] = {5.639,-1.685,1.803,-.3129e-1,-.6618e-2,.2352e-3}; | |
86 | Float_t Sigmap[4] = {.5668,-.2200e-1,.3657e-3,-.2201e-5}; | |
87 | Float_t Sigman[4] = {.4185,-.9798e-2,.1052e-3,-.4238e-6}; | |
88 | ||
89 | Float_t rnsp = SppvsB[0]+SppvsB[1]*b+SppvsB[2]*(b*b)+SppvsB[3]*(b*b*b)+ | |
90 | SppvsB[4]*(b*b*b*b)+SppvsB[5]*(b*b*b*b*b); | |
91 | Float_t rnsn = SpnvsB[0]+SpnvsB[1]*b+SpnvsB[2]*(b*b)+SpnvsB[3]*(b*b*b)+ | |
92 | SpnvsB[4]*(b*b*b*b)+SpnvsB[5]*(b*b*b*b*b); | |
93 | Float_t snsp = Sigmap[0]+Sigmap[1]*rnsp+Sigmap[2]*(rnsp*rnsp)+Sigmap[3]* | |
94 | (rnsp*rnsp*rnsp); | |
95 | Float_t snsn = Sigman[0]+Sigman[1]*rnsn+Sigman[2]*(rnsn*rnsn)+Sigman[3]* | |
96 | (rnsn*rnsn*rnsn); | |
97 | ||
98 | snsp = snsp*rnsp; | |
99 | snsn = snsn*rnsn; | |
100 | ||
101 | Float_t xgaup = gRandom->Gaus(0.0,1.0); | |
102 | snsp = snsp*xgaup; | |
103 | Float_t xgaun = gRandom->Gaus(0.0,1.0); | |
104 | snsn = snsn*xgaun; | |
105 | rnsp=rnsp+snsp; | |
106 | rnsn=rnsn+snsn; | |
107 | ||
108 | *NSpectp = Int_t(rnsp); | |
109 | *NSpectn = Int_t(rnsn); | |
110 | ||
111 | } |