Added AliGenHepMCEventHeader class and put it in STEERBase library. Updated AliGenRea...
[u/mrichter/AliRoot.git] / EVGEN / AliGenEMlib.cxx
index 5be02c9..0548554 100644 (file)
@@ -40,39 +40,96 @@ Int_t AliGenEMlib::fgSelectedPtParam=AliGenEMlib::kPizero7TeVpp;
 Int_t AliGenEMlib::fgSelectedCentrality=AliGenEMlib::kpp;
 Int_t AliGenEMlib::fgSelectedV2Systematic=AliGenEMlib::kNoV2Sys;
 
-Double_t AliGenEMlib::CrossOverLc(const double a, const double b, const double x){
+Double_t AliGenEMlib::CrossOverLc(double a, double b, double x){
   if(x<b-a/2) return 1.0;
   else if(x>b+a/2) return 0.0;
   else return cos(((x-b)/a+0.5)*TMath::Pi())/2+0.5;
 }
-Double_t AliGenEMlib::CrossOverRc(const double a, const double b, const double x){
+Double_t AliGenEMlib::CrossOverRc(double a, double b, double x){
   return 1-CrossOverLc(a,b,x);
 }
 
-const Double_t AliGenEMlib::fgkV2param[16][15] = {
-  // charged pion                                                                                                                        cent, based on
-  {  0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 }   // pp no V2
-  ,{ 6.792545e-02, 1.324227e+00, 2.071866e-02, 3.385157e+00, 5.183718e+00, 2.783082e-02, 5.000000e+00, 4.078812e+00, 5.444857e-02, 1.555569e-03, 1.099253e+00, 0, 1, 8.823293e-04, 1.310387e+00 }   // 0-5,  Francesco
-  ,{ 1.219655e-01, 1.237698e+00, 2.595132e-02, 3.157085e+00, 3.697252e+00, 4.576262e-02, 1.103851e+00, 5.102592e+00, 1.072299e-01, 2.338602e-03, 8.475128e-01, 0, 1, 1.631554e-03, 9.432292e-01 }   // 5-10, Francesco
-  ,{ 1.758981e-01, 1.253238e+00, 2.823197e-02, 3.455476e+00, 4.283204e+00, 6.053191e-02, 6.737159e-01, 4.420205e+00, 1.733567e-01, 3.292614e-03, 6.389341e-01, 0, 1, 3.234254e-03, 5.304084e-01 }   // 10-20,Francesco
-  ,{ 2.231108e-01, 1.316781e+00, 3.537289e-02, 3.228243e+00, 3.697226e+00, 7.131631e-02, 5.902192e-01, 4.938186e+00, 2.123481e-01, 4.775233e-03, 5.707533e-01, 0, 1, 3.906720e-03, 5.261247e-01 }   // 20-30,Francesco
-  ,{ 2.045341e-01, 1.324054e+00, 1.125395e-01, 2.179936e+00, 5.614603e+00, 8.946014e-02, 5.601938e-01, 2.959228e+00, 3.669442e-01, 5.318697e-03, 6.032741e-01, 0, 1, 4.231316e-03, 5.563485e-01 }   // 30-40,Francesco
-  ,{ 2.074778e-01, 1.894322e+00, 6.019224e-02, 2.414500e+00, 3.822721e+00, 9.891310e-02, 4.976126e-01, 1.256500e+00, 5.210237e-01, 5.687290e-03, 6.785355e-01, 0, 1, 3.954051e-03, 6.913428e-01 }   // 40-50,Francesco
-  ,{ 1.961099e-01, 2.122719e+00, 7.633636e-02, 1.731213e+00, 2.427536e+00, 1.026013e-01, 7.932261e-01, 4.082226e+00, 2.495494e-01, 5.677386e-03, 8.405024e-01, 0, 1, 2.961753e-03, 1.008892e+00 }   // 50-60,Francesco
-  ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 }   // 0-10
-  ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 }   // 20-40
-  ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 }   // 40-60
-  ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 }   // 60-80
-  ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 }   // 0-20
-  ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 }   // 0-40
-  ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 }   // 20-80
-  ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 }   // 40-80
+const Double_t AliGenEMlib::fgkV2param[kCentralities][16] = {
+  // charged pion                                                                                                                        cent, based on: https://twiki.cern.ch/twiki/bin/viewauth/ALICE/FlowPAGQM2012talkIdentified
+  {  0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000, 10 }  // pp no V2
+  ,{ 6.554571e-02, 1.436915e+00, 4.610598e-02, 2.554090e+00, 1.300948e+00, 2.970850e-02, 4.767877e+00, 4.228885e+00, 6.025959e-02, 1.570851e-03, 1.108941e+00, 0, 1, 1.715434e-03, 4.088070e-01, 25 }  // 0-5
+  ,{ 1.171348e-01, 1.333067e+00, 4.537086e-02, 3.046348e+00, 3.903416e+00, 4.407152e-02, 9.123846e-01, 4.834531e+00, 1.186227e-01, 2.259463e-03, 8.916458e-01, 0, 1, 2.300647e-03, 4.531172e-01, 25 }  // 5-10
+  ,{ 1.748434e-01, 1.285199e+00, 4.219881e-02, 4.018858e+00, 4.255082e+00, 7.955896e-03, 1.183264e-01,-9.329627e+00, 5.826570e-01, 3.368057e-03, 5.437668e-01, 0, 1, 3.178663e-03, 3.617552e-01, 25 }  // 10-20
+  ,{ 2.149526e-01, 1.408792e+00, 5.062101e-02, 3.206279e+00, 3.988517e+00, 3.724655e-02, 1.995791e-01,-1.571536e+01, 6.494227e+00, 4.957874e-03, 4.903140e-01, 0, 1, 4.214626e-03, 3.457922e-01, 25 }  // 20-30
+  ,{ 2.408942e-01, 1.477541e+00, 5.768983e-02, 3.333347e+00, 3.648508e+00,-2.044309e-02, 1.004145e-01,-2.386625e+01, 3.301913e+00, 5.666750e-03, 5.118686e-01, 0, 1, 4.626802e-03, 3.188974e-01, 25 }  // 30-40
+  ,{ 2.495109e-01, 1.543680e+00, 6.217835e-02, 3.518863e+00, 4.557145e+00, 6.014553e-02, 1.491814e-01,-5.443647e+00, 5.403300e-01, 6.217285e-03, 5.580218e-01, 0, 1, 4.620486e-03, 3.792879e-01, 25 }  // 40-50
+  ,{ 2.166399e-01, 1.931033e+00, 8.195390e-02, 2.226640e+00, 3.106649e+00, 1.058755e-01, 8.557791e-01, 4.006501e+00, 2.476449e-01, 6.720714e-03, 6.342966e-01, 0, 1, 4.449839e-03, 4.968750e-01, 25 }  // 50-60
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000, 10 }  // 0-10 
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000, 10 }  // 20-40
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000, 10 }  // 40-60
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000, 10 }  // 60-80
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000, 10 }  // 0-20 
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000, 10 }  // 0-40 
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000, 10 }  // 20-80
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000, 10 }  // 40-80
+};
+
+const Double_t AliGenEMlib::fgkRawPtOfV2Param[kCentralities][10] = {
+  {  0.0000000000, 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000 } // pp no V2
+  ,{ 2.181446e+08, 9.412925e-01, 1.158774e-01, 3.020303e+01, 6.790828e+00, 9.999996e+01, 2.616827e+00, 3.980492e+00, 1.225169e+07, 5.575243e+00 } // 0-5
+  ,{ 3.006215e+08, 9.511881e-01, 1.192788e-01, 2.981931e+01, 5.068175e+01, 9.999993e+01, 2.650635e+00, 4.073982e+00, 2.508045e+07, 5.621039e+00 } // 5-10
+  ,{ 1.643438e+09, 9.604242e-01, 1.218512e-01, 2.912684e+01, 1.164242e+00, 9.999709e+01, 2.662326e+00, 4.027795e+00, 7.020810e+07, 5.696860e+00 } // 10-20
+  ,{ 8.109985e+08, 9.421935e-01, 1.328020e-01, 2.655910e+01, 1.053677e+00, 9.999812e+01, 2.722949e+00, 3.964547e+00, 6.104096e+07, 5.694703e+00 } // 20-30
+  ,{ 5.219789e+08, 9.417339e-01, 1.417541e-01, 2.518080e+01, 7.430803e-02, 9.303295e+01, 2.780227e+00, 3.909570e+00, 4.723116e+07, 5.778375e+00 } // 30-40
+  ,{ 2.547159e+08, 9.481459e-01, 2.364858e-01, 1.689288e+01, 3.858883e+00, 6.352619e+00, 2.742270e+00, 3.855226e+00, 3.120535e+07, 5.878677e+00 } // 40-50
+  ,{ 9.396097e+07, 9.304491e-01, 3.244940e-01, 1.291696e+01, 2.854367e+00, 6.325908e+00, 2.828258e+00, 4.133699e+00, 1.302739e+07, 5.977896e+00 } // 50-60
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000 } // 0-10 
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000 } // 20-40
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000 } // 40-60
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000 } // 60-80
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000 } // 0-20 
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000 } // 0-40 
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000 } // 20-80
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000 } // 40-80
+};
+
+const Double_t AliGenEMlib::fgkPtParam[kCentralities][10] = {
+  {  0.0000000000, 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000 } // pp no V2
+  ,{ 7.641493e+01, 7.203468e-01, 3.651383e-01, 1.047542e+01, 3.494331e+00, 5.129019e+00, 3.081716e+00, 5.154525e+00, 3.065719e+01, 5.823718e+00 } // 0-5  
+  ,{ 1.704676e+02, 7.852682e-01, 4.177172e-01, 1.014652e+01, 3.324409e+00, 4.825894e+00, 2.889738e+00, 4.772249e+00, 3.502627e+01, 5.938773e+00 } // 5-10 
+  ,{ 1.823377e+02, 8.144309e-01, 4.291562e-01, 1.022767e+01, 3.585469e+00, 5.275078e+00, 3.144351e+00, 5.259097e+00, 2.675708e+01, 5.892506e+00 } // 10-20
+  ,{ 4.851407e+02, 9.341151e-01, 4.716673e-01, 1.058090e+01, 4.681218e+00, 7.261284e+00, 3.883227e+00, 6.638627e+00, 1.562806e+01, 5.772127e+00 } // 20-30
+  ,{ 3.157060e+01, 6.849451e-01, 4.868669e-01, 8.394558e+00, 3.539142e+00, 5.495280e+00, 4.102638e+00, 3.722991e+00, 1.638622e+01, 5.935963e+00 } // 30-40
+  ,{ 1.069397e+01, 5.816587e-01, 6.542961e-01, 6.472858e+00, 2.643870e+00, 3.929020e+00, 3.339224e+00, 2.410371e+00, 9.606748e+00, 6.116685e+00 } // 40-50
+  ,{ 1.857919e+01, 6.185989e-01, 5.878869e-01, 7.035064e+00, 2.892415e+00, 4.339383e+00, 3.549679e+00, 2.821061e+00, 1.529318e+01, 6.091388e+00 } // 50-60
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000 } // 0-10 
+  ,{ 1.271594e+02, 7.790165e-01, 5.793214e-01, 8.050008e+00, 3.211312e+00, 4.825258e+00, 3.840509e+00, 3.046231e+00, 2.172177e+01, 5.983496e+00 } // 20-40
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000 } // 40-60
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000 } // 60-80
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000 } // 0-20 
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000 } // 0-40 
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000 } // 20-80
+  ,{ 0.0000000000, 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000 } // 40-80
+};
+
+const Double_t AliGenEMlib::fgkThermPtParam[kCentralities][2] = {
+  {  0.0000000000, 0.0000000000 } // pp no V2
+  ,{ 0.0000000000, 0.0000000000 } // 0-5
+  ,{ 0.0000000000, 0.0000000000 } // 5-10
+  ,{ 3.335661e+01, 3.400585e+00 } // 10-20 //based on: https://aliceinfo.cern.ch/Notes/node/249
+  ,{ 0.0000000000, 0.0000000000 } // 20-30
+  ,{ 0.0000000000, 0.0000000000 } // 30-40
+  ,{ 0.0000000000, 0.0000000000 } // 40-50
+  ,{ 0.0000000000, 0.0000000000 } // 50-60
+  ,{ 3.648327e+02, 4.477749e+00 } // 0-10  //based on: https://aliceinfo.cern.ch/Notes/node/249
+  ,{ 1.696223e+00, 2.429660e+00 } // 20-40 //based on: https://twiki.cern.ch/twiki/pub/ALICE/ALICEDirectPhotonSpectrumPaper/directPbPb.pdf
+  ,{ 0.0000000000, 0.0000000000 } // 40-60
+  ,{ 0.0000000000, 0.0000000000 } // 60-80
+  ,{ 1.492160e+01, 2.805213e+00 } // 0-20  //based on: https://twiki.cern.ch/twiki/pub/ALICE/ALICEDirectPhotonSpectrumPaper/directPbPb.pdf
+  ,{ 4.215110e+01, 3.242719e+00 } // 0-40  //based on: https://aliceinfo.cern.ch/Figure/node/2866
+  ,{ 0.0000000000, 0.0000000000 } // 20-80
+  ,{ 0.0000000000, 0.0000000000 } // 40-80
 };
 
 // MASS   0=>PIZERO, 1=>ETA, 2=>RHO, 3=>OMEGA, 4=>ETAPRIME, 5=>PHI, 6=>JPSI
-const Double_t AliGenEMlib::fgkHM[7] = {0.13498, 0.54751, 0.7755, 0.78265, 0.95778, 1.01946, 3.0969};
+const Double_t AliGenEMlib::fgkHM[8] = {0.13498, 0.54751, 0.7755, 0.78265, 0.95778, 1.01946, 3.0969, 0.0};
 
-const Double_t AliGenEMlib::fgkMtFactor[2][7] = { 
+const Double_t AliGenEMlib::fgkMtFactor[2][8] = { 
   // {1.0, 0.5, 1.0, 0.9, 0.4, 0.23, 0.054},  // factor for pp from arXiv:1110.3929
   // {1.0, 0.55, 1.0, 0.9, 0.4, 0.25, 0.004}    // factor for PbPb from arXiv:1110.3929
   //{1., 0.48, 1.0, 0.9, 0.25, 0.4}, (old values)
@@ -81,10 +138,10 @@ const Double_t AliGenEMlib::fgkMtFactor[2][7] = {
   //https://aliceinfo.cern.ch/Figure/node/2634
   //https://aliceinfo.cern.ch/Figure/node/2788
   //https://aliceinfo.cern.ch/Figure/node/4403
-  //J/Psi PbPb from Comparison with Julian Books J/Psi -> e+e-, might be contradicting with https://aliceinfo.cern.ch/Figure/node/3457
+  //https://aliceinfo.cern.ch/Notes/node/87
   //best guess:
-  {1., 0.4, 1.0, 0.8, 0.4, 0.2, 0.0054}, //pp
-  {1., 0.4, 1.0, 0.8, 0.4, 0.25, 0.01}  //PbPb
+  {1., 0.48, 1.0, 0.9, 0.4, 0.25, 0., 0.}, //pp
+  {1., 0.48, 1.0, 0.9, 0.4, 0.25, 0., 0.}  //PbPb
 };
 
 //==========================================================================
@@ -98,13 +155,13 @@ const Double_t AliGenEMlib::fgkMtFactor[2][7] = {
 //                              General functions
 //
 //--------------------------------------------------------------------------
-Double_t AliGenEMlib::PtModifiedHagedornThermal(const Double_t pt,
-                                                const Double_t c,
-                                                const Double_t p0,
-                                                const Double_t p1,
-                                                const Double_t n,
-                                                const Double_t cT,
-                                                const Double_t T)
+Double_t AliGenEMlib::PtModifiedHagedornThermal(Double_t pt,
+                                                Double_t c,
+                                                Double_t p0,
+                                                Double_t p1,
+                                                Double_t n,
+                                                Double_t cT,
+                                                Double_t T)
 {
   // Modified Hagedorn Thermal fit to Picharged for PbPb:
   Double_t invYield;
@@ -115,12 +172,12 @@ Double_t AliGenEMlib::PtModifiedHagedornThermal(const Double_t pt,
 
 
 
-Double_t AliGenEMlib::PtModifiedHagedornExp(const Double_t pt,
-                                           const Double_t c,
-                                           const Double_t p1,
-                                           const Double_t p2,
-                                           const Double_t p0,
-                                           const Double_t n)
+Double_t AliGenEMlib::PtModifiedHagedornExp(Double_t pt,
+                                           Double_t c,
+                                           Double_t p1,
+                                           Double_t p2,
+                                           Double_t p0,
+                                           Double_t n)
 {
   // Modified Hagedorn exponentiel fit to Pizero for PbPb:
   Double_t invYield;
@@ -130,14 +187,14 @@ Double_t AliGenEMlib::PtModifiedHagedornExp(const Double_t pt,
 }
 
 
-Double_t AliGenEMlib::PtModifiedHagedornExp2(const Double_t pt,
-                                             const Double_t c,
-                                             const Double_t a,
-                                             const Double_t b,
-                                             const Double_t p0,
-                                             const Double_t p1,
-                                             const Double_t d,
-                                             const Double_t n)
+Double_t AliGenEMlib::PtModifiedHagedornExp2(Double_t pt,
+                                             Double_t c,
+                                             Double_t a,
+                                             Double_t b,
+                                             Double_t p0,
+                                             Double_t p1,
+                                             Double_t d,
+                                             Double_t n)
 {
   // Modified Hagedorn exponential fit to charged pions for pPb:
   Double_t invYield;
@@ -146,11 +203,11 @@ Double_t AliGenEMlib::PtModifiedHagedornExp2(const Double_t pt,
   return invYield*(2*TMath::Pi()*pt);
 }
 
-Double_t AliGenEMlib::PtTsallis(const Double_t pt,
-                                const Double_t m,
-                                const Double_t c,
-                                const Double_t T,
-                                const Double_t n)
+Double_t AliGenEMlib::PtTsallis(Double_t pt,
+                                Double_t m,
+                                Double_t c,
+                                Double_t T,
+                                Double_t n)
 {
   // Tsallis fit to Pizero for pp:
   Double_t mt;
@@ -162,7 +219,123 @@ Double_t AliGenEMlib::PtTsallis(const Double_t pt,
   return invYield*(2*TMath::Pi()*pt);
 }
 
+// Exponential
+Double_t AliGenEMlib::PtExponential(const Double_t *px, const Double_t *c){
+  const double &pt=px[0];
+  Double_t invYield = c[0]*exp(-pt*c[1]);
+  
+  return invYield*(2*TMath::Pi()*pt);
+}
+
+// Hagedorn with additional Powerlaw
+Double_t AliGenEMlib::PtModifiedHagedornPowerlaw(const Double_t *px, const Double_t *c){
+  const double &pt=px[0];
+  Double_t invYield = c[0]*pow(c[1]+pt*c[2],-c[3])*CrossOverLc(c[5],c[4],pt)+CrossOverRc(c[7],c[6],pt)*c[8]*pow(pt+0.001,-c[9]); //pt+0.001: prevent powerlaw from exploding for pt->0
+  
+  return invYield*(2*TMath::Pi()*pt+0.001); //pt+0.001: be sure to be > 0
+}
+
+// double powerlaw for J/Psi yield
+Double_t AliGenEMlib::PtDoublePowerlaw(const Double_t *px, const Double_t *c){
+  const double &pt=px[0];
+  Double_t yield = c[0]*pt*pow(1+pow(pt*c[1],2),-c[2]);
+  
+  return yield;
+}
+
+// integral over krollwada with S=1^2*(1-mee^2/mh^2)^3 from mee=0 up to mee=mh
+// approximation is perfect for mh>20MeV
+Double_t AliGenEMlib::IntegratedKrollWada(const Double_t *mh, const Double_t *){
+  if(*mh<0.002941) return 0;
+  return 2*log(*mh/0.000511/exp(1.75))/411.11/TMath::Pi();
+}
+
+//--------------------------------------------------------------------------
+//
+//                             DirectRealGamma
+//
+//--------------------------------------------------------------------------
+Double_t AliGenEMlib::PtPromptRealGamma( const Double_t *px, const Double_t */*dummy*/ )
+{
+  const static Double_t promptGammaPtParam[10] = { 8.715017e-02, 4.439243e-01, 1.011650e+00, 5.193789e+00, 2.194442e+01, 1.062124e+01, 2.469876e+01, 6.052479e-02, 5.611410e-02, 5.169743e+00 };
+  return PtModifiedHagedornPowerlaw(px,promptGammaPtParam)*GetTAA(fgSelectedCentrality);
+}
+
+Double_t AliGenEMlib::PtThermalRealGamma( const Double_t *px, const Double_t */*dummy*/ )
+{
+  return PtExponential(px,fgkThermPtParam[fgSelectedCentrality]);
+}
+
+Double_t AliGenEMlib::PtDirectRealGamma( const Double_t *px, const Double_t */*dummy*/ )
+{
+  return PtPromptRealGamma(px,px)+PtThermalRealGamma(px,px);
+}
+
+Int_t AliGenEMlib::IpDirectRealGamma(TRandom *)
+{
+  return 22;
+}
+
+Double_t AliGenEMlib::YDirectRealGamma( const Double_t *px, const Double_t */*dummy*/ )
+{
+  return YFlat(*px);
+}
+
+Double_t AliGenEMlib::V2DirectRealGamma( const Double_t *px, const Double_t */*dummy*/ )
+{
+  const static Double_t v2Param[3][16] = {
+    { 1.004245e-01, 1.057645e+00, 0.000000e+00, 2.836492e+00, 2.819767e+00, -6.231529e-02, 1.173054e+00, 2.836492e+00, 1.881590e-01, 1.183293e-02, 1.252249e+00, 0, 1, 4.876263e-03, 1.518526e+00, 4.5 } // 00-20, based on: https://aliceinfo.cern.ch/Notes/node/249
+    ,{ 1.619000e-01, 1.868201e+00, 6.983303e-15, 2.242170e+00, 4.484339e+00, -1.695734e-02, 2.301359e+00, 2.871469e+00, 1.619000e-01, 2.264320e-02, 1.028641e+00, 0, 1, 8.172203e-03, 1.271637e+00, 4.5 } // 20-40
+    ,{ 1.335000e-01, 1.076916e+00, 1.462605e-08, 2.785732e+00, 5.571464e+00, -2.356156e-02, 2.745437e+00, 2.785732e+00, 1.335000e-01, 1.571589e-02, 1.001131e+00, 0, 1, 5.179715e-03, 1.329344e+00, 4.5 } // 00-40
+  };
+  switch(fgSelectedCentrality){
+  case k0020: return V2Param(px,v2Param[0]); break;
+  case k2040: return V2Param(px,v2Param[1]); break;
+  case k0040: return V2Param(px,v2Param[2]); break;
+    // case k0010: return 0.43*V2Param(px,v2Param[1]); break;  //V2Pizero(0010)/V2Pizero(2040)=0.43 +-0.025
+    // case k1020: return 0.75*V2Param(px,v2Param[1]); break;  //V2Pizero(1020)/V2Pizero(2040)=0.75 +-0.04
+  case k0010: return 0.53*V2Param(px,v2Param[2]); break;  //V2Pizero(0010)/V2Pizero(0040)=0.53 +-0.03
+  case k1020: return 0.91*V2Param(px,v2Param[2]); break;  //V2Pizero(1020)/V2Pizero(0040)=0.91 +-0.04
+  }
+  return 0;
+}
+
+
+//--------------------------------------------------------------------------
+//
+//                             DirectVirtGamma
+//
+//--------------------------------------------------------------------------
+Double_t AliGenEMlib::PtPromptVirtGamma( const Double_t *px, const Double_t */*dummy*/ )
+{
+  return IntegratedKrollWada(px,px)*PtPromptRealGamma(px,px);
+}
+
+Double_t AliGenEMlib::PtThermalVirtGamma( const Double_t *px, const Double_t */*dummy*/ )
+{
+  return IntegratedKrollWada(px,px)*PtThermalRealGamma(px,px);
+}
+
+Double_t AliGenEMlib::PtDirectVirtGamma( const Double_t *px, const Double_t */*dummy*/ )
+{
+  return IntegratedKrollWada(px,px)*PtDirectRealGamma(px,px);
+}
+
+Int_t AliGenEMlib::IpDirectVirtGamma(TRandom *)
+{
+  return 220000;
+}
+
+Double_t AliGenEMlib::YDirectVirtGamma( const Double_t *px, const Double_t */*dummy*/ )
+{
+  return YFlat(*px);
+}
 
+Double_t AliGenEMlib::V2DirectVirtGamma( const Double_t *px, const Double_t */*dummy*/ )
+{
+  return V2DirectRealGamma(px,px);
+}
 
 //--------------------------------------------------------------------------
 //
@@ -177,6 +350,13 @@ Int_t AliGenEMlib::IpPizero(TRandom *)
 
 Double_t AliGenEMlib::PtPizero( const Double_t *px, const Double_t */*dummy*/ )
 {
+  // double pigammacorr=1; //misuse pion for direct gammas, tuned for 0040, iteration 0
+  // pigammacorr*=2.258900e-01*log(*px+0.001)+1.591291e+00;  //iteration 1
+  // pigammacorr*=6.601943e-03*log(*px+0.001)+9.593698e-01;  //iteration 2
+  // pigammacorr*=4.019933e-03*log(*px+0.001)+9.843412e-01;  //iteration 3
+  // pigammacorr*=-4.543991e-03*log(*px+0.001)+1.010886e+00; //iteration 4
+  // return pigammacorr*PtPromptRealGamma(px,px); //now the gammas from the pi->gg decay have the pt spectrum of prompt real gammas
+  
   // fit functions and corresponding parameter of Pizero pT for pp @ 2.76 TeV and @ 7 TeV and for PbPb @ 2.76 TeV 
 
   Double_t km=0.;
@@ -191,56 +371,115 @@ Double_t AliGenEMlib::PtPizero( const Double_t *px, const Double_t */*dummy*/ )
   Double_t kb=0.;
   Double_t kd=0.;
 
+  double n1,n2,n3,n4;
+  int oldCent;
+
   switch(fgSelectedPtParam|fgSelectedCentrality) {
+    // fit to pi charged, same data like in kPiOldChargedPbPb,
+    // but tested and compared against newest (2014) neutral pi measurement
+  case kPichargedPbPb|k0005:
+  case kPichargedPbPb|k0510:
+  case kPichargedPbPb|k1020:
+  case kPichargedPbPb|k2030:
+  case kPichargedPbPb|k3040:
+  case kPichargedPbPb|k4050:
+  case kPichargedPbPb|k5060:
+  case kPichargedPbPb|k2040:
+    return PtModifiedHagedornPowerlaw(px,fgkPtParam[fgSelectedCentrality]);
+    break;
+  case kPichargedPbPb|k0010:
+    n1=PtModifiedHagedornPowerlaw(px,fgkPtParam[k0005]);
+    n2=PtModifiedHagedornPowerlaw(px,fgkPtParam[k0510]);
+    return (n1+n2)/2;
+    break;
+  case kPichargedPbPb|k0020:
+    n1=PtModifiedHagedornPowerlaw(px,fgkPtParam[k0005]);
+    n2=PtModifiedHagedornPowerlaw(px,fgkPtParam[k0510]);
+    n3=PtModifiedHagedornPowerlaw(px,fgkPtParam[k1020]);
+    return (n1+n2+2*n3)/4;
+    break;
+  case kPichargedPbPb|k0040:
+    n1=PtModifiedHagedornPowerlaw(px,fgkPtParam[k0005]);
+    n2=PtModifiedHagedornPowerlaw(px,fgkPtParam[k0510]);
+    n3=PtModifiedHagedornPowerlaw(px,fgkPtParam[k1020]);
+    n4=PtModifiedHagedornPowerlaw(px,fgkPtParam[k2040]);
+    return (n1+n2+2*n3+4*n4)/8;
+    break;
+  case kPichargedPbPb|k4060:
+    n1=PtModifiedHagedornPowerlaw(px,fgkPtParam[k4050]);
+    n2=PtModifiedHagedornPowerlaw(px,fgkPtParam[k5060]);
+    return (n1+n2)/2;
+    break;
+
+
     // fit to pi charged v1
     // charged pion from ToF, unidentified hadrons scaled with pion from TPC
     // for Pb-Pb @ 2.76 TeV
-  case kPichargedPbPb|k0005:
+  case kPiOldChargedPbPb|k0005:
     kc=1347.5; kp0=0.9393; kp1=2.254; kn=11.294; kcT=0.002537; kT=2.414;
     return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
     break;
-  case kPichargedPbPb|k0510:
+  case kPiOldChargedPbPb|k0510:
     kc=1256.1; kp0=0.9545; kp1=2.248; kn=11.291; kcT=0.002662; kT=2.326;
     return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
     break;
-  case kPichargedPbPb|k2030:
+  case kPiOldChargedPbPb|k2030:
     kc=7421.6; kp0=1.2059; kp1=1.520; kn=10.220; kcT=0.002150; kT=2.196;
     return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
     break;
-  case kPichargedPbPb|k3040:
+  case kPiOldChargedPbPb|k3040:
     kc=1183.2; kp0=1.0478; kp1=1.623; kn=9.8073; kcT=0.00198333; kT=2.073;
     return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
     break;
     // the following is what went into the Pb-Pb preliminary approval (0-10%)
-  case kPichargedPbPb|k0010:
+  case kPiOldChargedPbPb|k0010:
     kc=1296.0; kp0=0.968; kp1=2.567; kn=12.27; kcT=0.004219; kT=2.207;
     return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
     break;
-  case kPichargedPbPb|k1020:
+  case kPiOldChargedPbPb|k1020:
     kc=986.0; kp0=0.9752; kp1=2.376; kn=11.62; kcT=0.003116; kT=2.213;
     return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
     break;
-  case kPichargedPbPb|k2040:
+  case kPiOldChargedPbPb|k2040:
     kc=17337.0; kp0=1.337; kp1=1.507; kn=10.629; kcT=0.00184; kT=2.234;
     return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
     break;
-  case kPichargedPbPb|k4050:
+  case kPiOldChargedPbPb|k4050:
     kc=6220.0; kp0=1.322; kp1=1.224; kn=9.378; kcT=0.000595; kT=2.383;
     return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
     break;
-  case kPichargedPbPb|k5060:
+  case kPiOldChargedPbPb|k5060:
     kc=2319.0; kp0=1.267; kp1=1.188; kn=9.044; kcT=0.000437; kT=2.276;
     return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
     break;
-  case kPichargedPbPb|k4060:
+  case kPiOldChargedPbPb|k4060:
     kc=4724.0; kp0=1.319; kp1=1.195; kn=9.255; kcT=0.000511; kT=2.344;
     return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
     break;
-  case kPichargedPbPb|k6080:
+  case kPiOldChargedPbPb|k6080:
     kc=2842.0; kp0=1.465; kp1=0.8324; kn=8.167; kcT=0.0001049; kT=2.29;
     return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
     break;
-   
+  case kPiOldChargedPbPb|k0020:
+    oldCent=fgSelectedCentrality;
+    fgSelectedCentrality=k0010;
+    n1=PtPizero(px,px);
+    fgSelectedCentrality=k1020;
+    n2=PtPizero(px,px);
+    fgSelectedCentrality=oldCent;
+    return (n1+n2)/2;
+    break;
+  case kPiOldChargedPbPb|k0040:
+    oldCent=fgSelectedCentrality;
+    fgSelectedCentrality=k0010;
+    n1=PtPizero(px,px);
+    fgSelectedCentrality=k1020;
+    n2=PtPizero(px,px);
+    fgSelectedCentrality=k2040;
+    n3=PtPizero(px,px);
+    fgSelectedCentrality=oldCent;
+    return (n1+n2+2*n3)/4;
+    break;
 
     // fit to pizero from conversion analysis
     // for PbPb @ 2.76 TeV
@@ -335,39 +574,51 @@ Double_t AliGenEMlib::YPizero( const Double_t *py, const Double_t */*dummy*/ )
 
 Double_t AliGenEMlib::V2Pizero( const Double_t *px, const Double_t */*dummy*/ )
 {
-  int selectedCentrality;
-  double n1,n2,v1,v2;
+  double n1,n2,n3,n4,n5;
+  double v1,v2,v3,v4,v5;
   switch(fgSelectedCentrality) {
   case k0010:
-    selectedCentrality=fgSelectedCentrality;
-    fgSelectedCentrality=k0005;
-    n1=PtPizero(px,(Double_t*) 0);
-    v1=V2Param(px,fgkV2param[fgSelectedCentrality]);
-    fgSelectedCentrality=k0510;
-    n2=PtPizero(px,(Double_t*) 0);
-    v2=V2Param(px,fgkV2param[fgSelectedCentrality]);
-    fgSelectedCentrality=selectedCentrality;
+    n1=PtModifiedHagedornPowerlaw(px,fgkRawPtOfV2Param[k0005]);
+    v1=V2Param(px,fgkV2param[k0005]);
+    n2=PtModifiedHagedornPowerlaw(px,fgkRawPtOfV2Param[k0510]);
+    v2=V2Param(px,fgkV2param[k0510]);
     return (n1*v1+n2*v2)/(n1+n2);
     break;
+  case k0020:
+    n1=PtModifiedHagedornPowerlaw(px,fgkRawPtOfV2Param[k0005]);
+    v1=V2Param(px,fgkV2param[k0005]);
+    n2=PtModifiedHagedornPowerlaw(px,fgkRawPtOfV2Param[k0510]);
+    v2=V2Param(px,fgkV2param[k0510]);
+    n3=PtModifiedHagedornPowerlaw(px,fgkRawPtOfV2Param[k1020]);
+    v3=V2Param(px,fgkV2param[k1020]);
+    return (n1*v1+n2*v2+2*n3*v3)/(n1+n2+2*n3);
+    break;
   case k2040:
-    selectedCentrality=fgSelectedCentrality;
-    fgSelectedCentrality=k2030;
-    n1=PtPizero(px,(Double_t*) 0);
-    v1=V2Param(px,fgkV2param[fgSelectedCentrality]);
-    fgSelectedCentrality=k3040;
-    n2=PtPizero(px,(Double_t*) 0);
-    v2=V2Param(px,fgkV2param[fgSelectedCentrality]);
-    fgSelectedCentrality=selectedCentrality;
+    n1=PtModifiedHagedornPowerlaw(px,fgkRawPtOfV2Param[k2030]);
+    v1=V2Param(px,fgkV2param[k2030]);
+    n2=PtModifiedHagedornPowerlaw(px,fgkRawPtOfV2Param[k3040]);
+    v2=V2Param(px,fgkV2param[k3040]);
     return (n1*v1+n2*v2)/(n1+n2);
     break;
+  case k0040:
+    n1=PtModifiedHagedornPowerlaw(px,fgkRawPtOfV2Param[k0005]);
+    v1=V2Param(px,fgkV2param[k0005]);
+    n2=PtModifiedHagedornPowerlaw(px,fgkRawPtOfV2Param[k0510]);
+    v2=V2Param(px,fgkV2param[k0510]);
+    n3=PtModifiedHagedornPowerlaw(px,fgkRawPtOfV2Param[k1020]);
+    v3=V2Param(px,fgkV2param[k1020]);
+    n4=PtModifiedHagedornPowerlaw(px,fgkRawPtOfV2Param[k2030]);
+    v4=V2Param(px,fgkV2param[k2030]);
+    n5=PtModifiedHagedornPowerlaw(px,fgkRawPtOfV2Param[k3040]);
+    v5=V2Param(px,fgkV2param[k3040]);
+    return (n1*v1+n2*v2+2*n3*v3+2*n4*v4+2*n5*v5)/(n1+n2+2*n3+2*n4+2*n5);
+    break;
 
   default:
     return V2Param(px,fgkV2param[fgSelectedCentrality]);
   }
 }
 
-
-
 //--------------------------------------------------------------------------
 //
 //                              Eta
@@ -561,7 +812,18 @@ Int_t AliGenEMlib::IpJpsi(TRandom *)
 Double_t AliGenEMlib::PtJpsi( const Double_t *px, const Double_t */*dummy*/ )
 {
   // Jpsi pT
-  return MtScal(*px,6);
+  // based on: //https://aliceinfo.cern.ch/Notes/node/242, https://aliceinfo.cern.ch/Figure/node/3457, www.sciencedirect.com/science/article/pii/S0370269312011446
+  const static Double_t jpsiPtParam[2][3] = {
+    {  9.686337e-03, 2.629441e-01, 4.552044e+00 }
+    ,{ 3.403549e-03, 2.897061e-01, 3.644278e+00 }
+  };
+  const double pt=px[0]*2.28/2.613;
+  switch(fgSelectedCentrality) {
+  case k0020: return 2.405*PtDoublePowerlaw(&pt,jpsiPtParam[0]); break;
+  case k2040: return 2.405*PtDoublePowerlaw(&pt,jpsiPtParam[1]); break;
+  case k0040: return 0.5*2.405*(PtDoublePowerlaw(&pt,jpsiPtParam[0])+PtDoublePowerlaw(&pt,jpsiPtParam[1])); break;
+  }
+  return 0;
 }
 
 Double_t AliGenEMlib::YJpsi( const Double_t *py, const Double_t */*dummy*/ )
@@ -571,18 +833,15 @@ Double_t AliGenEMlib::YJpsi( const Double_t *py, const Double_t */*dummy*/ )
 
 Double_t AliGenEMlib::V2Jpsi( const Double_t *px, const Double_t */*dummy*/ )
 {
-  const int oldSys=fgSelectedV2Systematic;
-  fgSelectedV2Systematic=kNoV2Sys;
-  double ret=0;
-
-  switch(oldSys){
-  case kLoV2Sys: ret=0; break;
-  case kNoV2Sys: ret=KEtScal(*px,6)/2; break;
-  case kUpV2Sys: ret=KEtScal(*px,6); break;
+  const static Double_t v2Param[16] = { 1.156000e-01, 8.936854e-01, 0.000000e+00, 4.000000e+00, 6.222375e+00, -1.600314e-01, 8.766676e-01, 7.824143e+00, 1.156000e-01, 3.484503e-02, 4.413685e-01, 0, 1, 3.484503e-02, 4.413685e-01, 7.2 };
+  switch(fgSelectedCentrality){
+  case k2040: return V2Param(px,v2Param); break;
+  case k0010: return 0.43*V2Param(px,v2Param); break;  //V2Pizero(0010)/V2Pizero(2040)=0.43 +-0.025
+  case k1020: return 0.75*V2Param(px,v2Param); break;  //V2Pizero(1020)/V2Pizero(2040)=0.75 +-0.04
+  case k0020: return 0.66*V2Param(px,v2Param); break;  //V2Pizero(0020)/V2Pizero(2040)=0.66 +-0.035
+  case k0040: return 0.82*V2Param(px,v2Param); break;  //V2Pizero(0040)/V2Pizero(2040)=0.82 +-0.05
   }
-
-  fgSelectedV2Systematic=oldSys;
-  return ret;
+  return 0;
 }
 
 Double_t AliGenEMlib::YFlat(Double_t /*y*/)
@@ -623,14 +882,11 @@ Double_t AliGenEMlib::MtScal(Double_t pt, Int_t np)
   return norm*(pt/scaledPt)*scaledYield;
 }
 
-Double_t AliGenEMlib::KEtScal(const Double_t pt, Int_t np)
+Double_t AliGenEMlib::KEtScal(Double_t pt, Int_t np)
 {
   const int nq=2; //number of quarks for particle np, here always 2
-  Double_t scaledPt = sqrt(pow(2/nq*(sqrt(pt*pt+fgkHM[np]*fgkHM[np])-fgkHM[np])+fgkHM[0],2)-fgkHM[0]*fgkHM[0]);
-  double val=V2Pizero(&scaledPt, (Double_t*) 0);
-  static const double syserr[12]={0., 0.09, 0.07, 0.06, 0.04, 0.04, 0.04, 0.05, 0., 0., 0., 0.}; //based on pi vs kaon
-  double sys=fgSelectedV2Systematic*min(fgkV2param[fgSelectedCentrality][0],fgkV2param[fgSelectedCentrality][8])*syserr[fgSelectedCentrality];
-  return TMath::Max(val+sys,0.0);
+  Double_t scaledPt = sqrt(pow(2.0/nq*(sqrt(pt*pt+fgkHM[np]*fgkHM[np])-fgkHM[np])+fgkHM[0],2)-fgkHM[0]*fgkHM[0]);
+  return V2Pizero(&scaledPt, (Double_t*) 0);
 }
 
 Double_t AliGenEMlib::V2Param(const Double_t *px, const Double_t *par)
@@ -639,18 +895,40 @@ Double_t AliGenEMlib::V2Param(const Double_t *px, const Double_t *par)
 
   const double &pt=px[0];
   double val=CrossOverLc(par[4],par[3],pt)*(2*par[0]/(1+TMath::Exp(par[1]*(par[2]-pt)))-par[0])+CrossOverRc(par[4],par[3],pt)*((par[8]-par[5])/(1+TMath::Exp(par[6]*(pt-par[7])))+par[5]);
-  double sys=fgSelectedV2Systematic*par[11+fgSelectedV2Systematic*2]*pow(pt,par[12+fgSelectedV2Systematic*2]);
-  return TMath::Max(val+sys,0.0);
+  double sys=0;
+  if(fgSelectedV2Systematic){
+    double syspt=pt>par[15]?par[15]:pt;
+    sys=fgSelectedV2Systematic*par[11+fgSelectedV2Systematic*2]*pow(syspt,par[12+fgSelectedV2Systematic*2]);
+  }
+  return std::max(val+sys,0.0);
 }
 
 Double_t AliGenEMlib::V2Flat(const Double_t */*px*/, const Double_t */*param*/)
 {
   // Flat v2
 
-  return 1.0;
+  return 0.0;
 }
 
-
+Double_t AliGenEMlib::GetTAA(Int_t cent){
+  const static Double_t taa[16] = { 1.0,    // pp
+                                   26.32,  // 0-5
+                                   20.56,  // 5-10
+                                   14.39,  // 10-20
+                                   8.70,   // 20-30
+                                   5.001,  // 30-40
+                                   2.675,  // 40-50
+                                   1.317,  // 50-60
+                                   23.44,  // 0-10
+                                   6.85,   // 20-40
+                                   1.996,  // 40-60
+                                   0.4174, // 60-80
+                                   18.91,  // 0-20
+                                   12.88,  // 0-40
+                                   3.088,  // 20-80
+                                   1.207}; // 40-80
+  return taa[cent];  
+}
 
 //==========================================================================
 //
@@ -670,6 +948,12 @@ GenFunc AliGenEMlib::GetPt(Int_t param, const char * tname) const
 
    switch (param) 
     {
+     case kDirectRealGamma:
+       func=PtDirectRealGamma;
+      break;
+     case kDirectVirtGamma:
+       func=PtDirectVirtGamma;
+      break;
     case kPizero:
       func=PtPizero;
       break;
@@ -707,6 +991,12 @@ GenFunc AliGenEMlib::GetY(Int_t param, const char * tname) const
 
    switch (param) 
     {
+    case kDirectRealGamma:
+      func=YDirectRealGamma;
+      break;
+    case kDirectVirtGamma:
+      func=YDirectVirtGamma;
+      break;
     case kPizero:
          func=YPizero;
          break;
@@ -744,6 +1034,12 @@ GenFuncIp AliGenEMlib::GetIp(Int_t param, const char * tname) const
 
    switch (param) 
     {
+    case kDirectRealGamma:
+      func=IpDirectRealGamma;
+      break;
+    case kDirectVirtGamma:
+      func=IpDirectVirtGamma;
+      break;
     case kPizero:
          func=IpPizero;
          break;
@@ -781,6 +1077,12 @@ GenFunc AliGenEMlib::GetV2(Int_t param, const char * tname) const
 
   switch (param) 
     {
+    case kDirectRealGamma:
+      func=V2DirectRealGamma;
+      break;
+    case kDirectVirtGamma:
+      func=V2DirectVirtGamma;
+      break;
     case kPizero:
       func=V2Pizero;
       break;