STARLIGHT update (missing modifications in 65f9e959b845e683a0d24545a6372e7f985b7853)

[u/mrichter/AliRoot.git] / STARLIGHT / starlight / TStarLight / testsl.C

diff --git a/STARLIGHT/starlight/TStarLight/testsl.C b/STARLIGHT/starlight/TStarLight/testsl.C
index eda0a5f683331f5573bc6e208d9b5e3bbcfd5d76..ce4b1af504a25a2830894b49dd7c9b8b041172ce 100644 (file)
--- a/STARLIGHT/starlight/TStarLight/testsl.C
+++ b/STARLIGHT/starlight/TStarLight/testsl.C
@@ -1,19 +1,21 @@
 // -*- C++ -*-
 // $Id$
 
+// creates histogram of generated rho0 kinematical data
+// writes a text file "rho0.txt" which should be identical to ../../test/rho0.txt
+
 void testsl() {
-  gSystem->AddDynamicPath("../../../objdir/lib/tgt_linux/");
   gSystem->Load("libStarLight");
   gSystem->Load("libAliStarLight.so");
 
   TStarLight* sl = new TStarLight("starlight generator", "title", "");
 
-  sl->SetParameter("BEAM_1_Z     =   82    #Z of projectile");
-  sl->SetParameter("BEAM_1_A     =  208    #A of projectile");
+  sl->SetParameter("BEAM_1_Z     =    1    #Z of projectile");
+  sl->SetParameter("BEAM_1_A     =    1    #A of projectile");
   sl->SetParameter("BEAM_2_Z     =   82    #Z of target");
   sl->SetParameter("BEAM_2_A     =  208    #A of target");
-  sl->SetParameter("BEAM_1_GAMMA = 1470    #Gamma of the colliding ions");
-  sl->SetParameter("BEAM_2_GAMMA = 1470    #Gamma of the colliding ions");
+  sl->SetParameter("BEAM_1_GAMMA = 4264.4  #Gamma of the colliding ions");
+  sl->SetParameter("BEAM_2_GAMMA = 1682.4  #Gamma of the colliding ions");
   sl->SetParameter("W_MAX        =   15    #Max value of w");
   sl->SetParameter("W_MIN        =    -1   #Min value of w");
   sl->SetParameter("W_N_BINS     =   40    #Bins i w");
@@ -26,32 +28,58 @@ void testsl() {
   sl->SetParameter("ETA_MIN      =  -10    #Minimum pseudorapidity");
   sl->SetParameter("ETA_MAX      =   10    #Maximum pseudorapidity");
   sl->SetParameter("PROD_MODE    =    2    #gg or gP switch (1 = 2-photon, 2 = coherent vector meson (narrow), 3 = coherent vector meson (wide), # 4 = incoherent vector meson, 5 = A+A DPMJet single, 6 = A+A DPMJet double, 7 = p+A DPMJet single, 8 = p+A Pythia single )");
-  // is N_EVENTS valid
-  sl->SetParameter("N_EVENTS     = 1000    #Number of events");
   sl->SetParameter("PROD_PID     =  113    #Channel of interest (not relevant for photonuclear processes)");
   sl->SetParameter("RND_SEED     = 34533   #Random number seed");
-  sl->SetParameter("OUTPUT_FORMAT =   2    #Form of the output");
   sl->SetParameter("BREAKUP_MODE  =   5    #Controls the nuclear breakup");
   sl->SetParameter("INTERFERENCE  =   0    #Interference (0 = off, 1 = on)");
   sl->SetParameter("IF_STRENGTH   =   1.   #% of intefernce (0.0 - 0.1)");
   sl->SetParameter("COHERENT      =   1    #Coherent=1,Incoherent=0");
   sl->SetParameter("INCO_FACTOR   =   1.   #percentage of incoherence");
-  sl->SetParameter("BFORD         =   9.5  #");
   sl->SetParameter("INT_PT_MAX    =   0.24 #Maximum pt considered, when interference is turned on");
   sl->SetParameter("INT_PT_N_BINS = 120    #Number of pt bins when interference is turned on");
 
   sl->InitStarLight();
   sl->PrintInputs(std::cout);
-
-  sl->GetParameter("INT_PT_MAX");
-
   TClonesArray tca("TParticle", 1000);
 
-  for (Int_t counter(0); counter<1000; ++counter) {
-    Printf("--------------------------------------------------------------------------------");
+  TLorentzVector v[2], vSum;
+  TH1* hM  = new TH1D("hM",  "STARLIGHT;M#(){#pi^{+}#pi^{-}}",     100, 0., 2.);
+  TH1* hPt = new TH1D("hPt", "STARLIGHT;P_{T}#(){#pi^{+}#pi^{-}}", 100, 0., 1.);
+  TH1* hY  = new TH1D("hY",  "STARLIGHT;Y#(){#pi^{+}#pi^{-}}",     160,-8., 8.);
+
+  std::ofstream ofs("rho0.txt");
+  TParticle *p;
+  for (Int_t counter(0); counter<100;) {
     sl->GenerateEvent();    
+    sl->BoostEvent();    
     sl->ImportParticles(&tca, "ALL");
-    for (Int_t i=0; i<tca.GetEntries(); ++i)
-      tca.At(i)->Print();
+    Bool_t isOK = kTRUE;
+    for (Int_t i=0; i<tca.GetEntries() && isOK; ++i) {
+      p = (TParticle*)tca.At(i);
+      p->Momentum(v[i]);
+//       isOK = TMath::Abs(v[i].Rapidity()) <= 0.9;
+    }
+    if (!isOK) continue;
+    Printf("counter, %d", counter, isOK);
+    ++counter;
+    vSum = v[0] + v[1];
+    ofs << std::fixed << std::setprecision(4)
+       << vSum.M() << " " << vSum.Perp() << " " << vSum.Rapidity() << " "
+       << v[0].Eta() << " " << v[0].Px() << " " << v[0].Py() << " " << v[0].Pz() << " "
+       << v[1].Eta() << " " << v[1].Px() << " " << v[1].Py() << " " << v[1].Pz()
+       << std::endl;
+    hM->Fill(vSum.M());
+    hPt->Fill(vSum.Perp());
+    hY->Fill(vSum.Rapidity());
   }
+  TFile::Open("sl.root", "RECREATE");
+  sl->Write();
+  gFile->Write();
+
+  hM->Draw();
+  c1->SaveAs("SL.pdf(");
+  hPt->Draw();
+  c1->SaveAs("SL.pdf");
+  hY->Draw();
+  c1->SaveAs("SL.pdf)");
 }