Hydjet2 update

GeneratorInterface/Hydjet2Interface/interface/DatabasePDG.h

@@ -27,7 +27,8 @@
 #include "FWCore/Utilities/interface/InputTag.h"
 #include "FWCore/ServiceRegistry/interface/Service.h"
 #include "CommonTools/UtilAlgos/interface/TFileService.h"
-#include "DataFormats/HepMCCandidate/interface/GenParticle.h"
+#include "DataFormats/HepMCCandidate/interface/GenParticle.h"
+
 const int kMaxParticles = 1000;
 
 class DatabasePDG {

GeneratorInterface/Hydjet2Interface/interface/DecayChannel.h

@@ -14,6 +14,7 @@
 #define DECAY_CHANNEL
 
 #include "Rtypes.h"
+#include "FWCore/MessageLogger/interface/MessageLogger.h"
 
 const int kMaxDaughters = 3;
 const int kNonsensePDG = 1000000000;

GeneratorInterface/Hydjet2Interface/interface/RandArrayFunction.h

@@ -2,7 +2,7 @@
 #define RANDARRAYFUNCTION_INCLUDED
 
 #include <vector>
-
+#include "FWCore/MessageLogger/interface/MessageLogger.h"
 #include "CLHEP/Random/RandomEngine.h"
 #include "CLHEP/Random/RandFlat.h"
 extern CLHEP::HepRandomEngine* hjRandomEngine;

GeneratorInterface/Hydjet2Interface/python/hydjet2DefaultParameters_cff.py

@@ -14,6 +14,17 @@
 	fR 	= cms.double(13.45), 	# Maximal transverse radius at thermal freeze-out for central collisions, fR [fm]
 )
 
+collisionParameters5100GeV = cms.PSet(
+	fAw = cms.double(208.0), ## beam/target atomic number
+	fSqrtS = cms.double(5100.0),
+	fUmax 	= cms.double(1.35), 	# Maximal transverse flow rapidity at thermal freeze-out for central collisions, fUmax
+	fPtmin 	= cms.double(10.), 	# Minimal pt of parton-parton scattering in PYTHIA event, fPtmin [GeV/c] 
+	fT0 = cms.double(1.1), 		# Initial QGP temperature for central Pb+Pb collisions in mid-rapidity, fT0 [GeV]; allowed range [0.2,2.0]GeV;
+
+	### Volume parameters at thermal freeze-out ###
+	fTau 	= cms.double(13.2), 	# Proper time proper at thermal freeze-out for central collisions, fTau [fm/c]
+	fR 	= cms.double(13.9), 	# Maximal transverse radius at thermal freeze-out for central collisions, fR [fm]
+)
 
 qgpParameters = cms.PSet(
 	fTau0 	= cms.double(0.1), 	# Proper QGP formation time in fm/c, fTau0 (0.01<fTau0<10)
@@ -60,12 +71,12 @@
 
 	fEtaType 	= cms.double(1.), 	# Flag to choose longitudinal flow rapidity distribution, fEtaType (=0 uniform, >0 Gaussian with the dispersion Ylmax)
 	fIshad 	= cms.int32(1), 	# Flag to switch on/off nuclear shadowing, fIshad (0 shadowing off, 1 shadowing on)	
-	fPyhist 	= cms.int32(1), 	# Flag to suppress the output of particle history from PYTHIA, fPyhist (=1 only final state particles; =0 full particle history from PYTHIA)
+	fPyhist 	= cms.int32(0), 	# Flag to suppress the output of particle history from PYTHIA, fPyhist (=1 only final state particles; =0 full particle history from PYTHIA)
 	fIenglu 	= cms.int32(0), 	# Flag to fix type of partonic energy loss, fIenglu (0 radiative and collisional loss, 1 radiative loss only, 2 collisional loss only)
 	fIanglu 	= cms.int32(1), 	# Flag to fix type of angular distribution of in-medium emitted gluons, fIanglu (0 small-angular, 1 wide-angular, 2 collinear).
 	embeddingMode = cms.bool(False),
 	rotateEventPlane = cms.bool(True)
-	#fSeed	= cms.int32(0) 	# Parameter to set random number seed, fSeed (=0 the current time is used, >0 the value fSeed is used) 
+
 )
 
 PythiaDefaultBlock = cms.PSet(

GeneratorInterface/Hydjet2Interface/python/hydjet2Default_cfi.py

@@ -3,7 +3,7 @@
 from GeneratorInterface.Hydjet2Interface.hydjet2DefaultParameters_cff import *
 
 generator = cms.EDFilter("Hydjet2GeneratorFilter",
-	collisionParameters2760GeV,
+	collisionParameters5100GeV,
 	qgpParameters,
 	hydjet2Parameters,
 	fNhsel 	= cms.int32(2), 	# Flag to include jet (J)/jet quenching (JQ) and hydro (H) state production, fNhsel (0 H on & J off, 1 H/J on & JQ off, 2 H/J/HQ on, 3 J on & H/JQ off, 4 H off & J/JQ on)

GeneratorInterface/Hydjet2Interface/src/DatabasePDG.cc

@@ -70,29 +70,28 @@ bool DatabasePDG::LoadData() {
 bool DatabasePDG::LoadParticles() {
   ifstream particleFile;
   particleFile.open(fParticleFilename);
-  if(!particleFile) {
-    cout << "ERROR in DatabasePDG::LoadParticles() : The ASCII file containing the PDG particle list (\""
-         << fParticleFilename << "\") was not found !! Aborting..." << endl;
+  if(!particleFile) {
+    edm::LogError("DatabasePDG")<< "The ASCII file containing the PDG particle list " << fParticleFilename << " was not found";
     return kFALSE;
   }
 
   char name[9];
   double mass, width, spin, isospin, isospinZ, q, s, aq, as, c, ac;
   int pdg;
-  int goodStatusParticles = 0;
+  int goodStatusParticles = 0;
 
-  cout << "Info in DatabasePDG::LoadParticles() : Start loading particles with the following criteria:" << endl
-       << "       Use particles containing charm quarks (1-yes;0-no) : " << fUseCharmParticles << endl
+  edm::LogInfo("DatabasePDG")<< "Start loading particles with the following criteria:" << endl
+       << "       Use particles containing charm quarks (1:yes;0:no) : " << fUseCharmParticles << endl
        << "       Mass range                                         : (" << fMinimumMass << "; " << fMaximumMass << ")" << endl
-       << "       Width range                                        : (" << fMinimumWidth << "; " << fMaximumWidth << ")" << endl;
+       << "       Width range                                        : (" << fMinimumWidth << "; " << fMaximumWidth << ")";
 
   particleFile.exceptions(ios::failbit);
   while(!particleFile.eof()) {
     try {
       particleFile >> name >> mass >> width >> spin >> isospin >> isospinZ >> q >> s >> aq >> as >> c >> ac >> pdg;
     }
     catch (ios::failure const &problem) {
-      cout << problem.what() << endl;
+      LogDebug("DatabasePDG")<<" ios:failure in particle file "<< problem.what();
       break;
     }
 
@@ -130,22 +129,21 @@ bool DatabasePDG::LoadParticles() {
     fNParticles++;
   }
   particleFile.close();
-  if(fNParticles==0) {
-    cout << "Warning in DatabasePDG::LoadParticles(): No particles were found in the file specified!!" << endl;
+  if(fNParticles==0) {
+
+    LogWarning("DatabasePDG")<<" No particles were found in the file specified!!";
     return kFALSE;
   }
   SortParticles();
-  cout << "Info in DatabasePDG::LoadParticles(): Particle definitions found = " << fNParticles << endl
-       << "                                      Good status particles      = " << goodStatusParticles << endl;
+  edm::LogInfo("DatabasePDG")<< " Particle definitions found: " << fNParticles << ". Good status particles: " << goodStatusParticles;
   return kTRUE;
 }
 
 bool DatabasePDG::LoadDecays() {
   ifstream decayFile;
   decayFile.open(fDecayFilename);
-  if(!decayFile) {
-    cout << "ERROR in DatabasePDG::LoadDecays() : The ASCII file containing the decays list (\""
-         << fDecayFilename << "\") was not found !! Aborting..." << endl;
+  if(!decayFile) {
+    edm::LogError("DatabasePDG")<< "The ASCII file containing the decays list " << fDecayFilename << " was not found";
     return kFALSE;
   }
 
@@ -163,8 +161,8 @@ bool DatabasePDG::LoadDecays() {
         decayFile >> daughter_pdg[i];
       decayFile >> branching;
     }
-    catch (ios::failure const &problem) {
-      cout << problem.what() << endl;
+    catch (ios::failure const &problem) {
+      LogDebug("DatabasePDG")<<" ios:failure in decay file "<< problem.what();
       break;
     }
     if((mother_pdg!=0) && (daughter_pdg[0]!=0) && (branching>=0)) {
@@ -174,20 +172,18 @@ bool DatabasePDG::LoadDecays() {
           nDaughters++;
       ParticlePDG* particle = GetPDGParticle(mother_pdg);
       if(!particle) {
-	cout << "DatabasePDG::LoadDecays(): WARNING!!! Mother particle PDG ("
-	     << mother_pdg << ") not found in the particle definition list: " << mother_pdg << "--> ";
+	LogWarning("DatabasePDG")<<" Mother particle PDG (" << mother_pdg 
+		<< ") not found in the particle definition list:"<< mother_pdg << " >>> ";
 	for(int kk=0; kk<nDaughters; kk++) 
-	  cout << daughter_pdg[kk] << "  ";
-	cout << endl;
+	  LogWarning("DatabasePDG")<< daughter_pdg[kk] << "  ";
 	return kFALSE;
       }
       for(int kk=0; kk<nDaughters; kk++) {
 	if(!GetPDGParticle(daughter_pdg[kk])) {
-	  cout << "DatabasePDG::LoadDecays(): WARNING!!! Daughter particle PDG ("
-	       << daughter_pdg[kk] << ") not found in the particle definition list: " << mother_pdg << "--> ";
+	  LogWarning("DatabasePDG")<<"Daughter particle PDG (" << daughter_pdg[kk] 
+		<< ") not found in the particle definition list: " << mother_pdg << ">>> ";
 	  for(int kkk=0; kkk<nDaughters; kkk++) 
-	    cout << daughter_pdg[kkk] << "  ";
-	  cout << endl;
+	    LogWarning("DatabasePDG")<< daughter_pdg[kkk] << "  ";
 	}
       }
       DecayChannel decay(mother_pdg, branching, nDaughters, daughter_pdg);
@@ -199,25 +195,25 @@ bool DatabasePDG::LoadDecays() {
   for(int i=0; i<fNParticles; i++) {
     nDecayChannels += fParticles[i]->GetNDecayChannels();
   }
-  cout << "Info in DatabasePDG::LoadDecays(): Number of decays found in the database is " << nDecayChannels << endl;
+  edm::LogInfo("DatabasePDG")<< "Number of decays found in the database is " << nDecayChannels;
   return kTRUE;
 }
 
 ParticlePDG* DatabasePDG::GetPDGParticleByIndex(int index) {
   if(index<0 || index>fNParticles) {
-    cout << "Warning in DatabasePDG::GetPDGParticleByIndex(int): Particle index is negative or too big !!" << endl
-         << " It must be inside this range: [0, " << fNParticles-1 << "]" << endl
-         << " Returning null pointer!!" << endl;
+    edm::LogWarning("DatabasePDG")<< "Particle index is negative or too big !!" << endl
+         << " It must be inside this range: (0, " << fNParticles-1 << ")" << endl
+         << " Returning null pointer!!";
     return 0x0;
   }
   return fParticles[index];
 }
 
 bool DatabasePDG::GetPDGParticleStatusByIndex(int index) {
   if(index<0 || index>fNParticles) {
-    cout << "Warning in DatabasePDG::GetPDGParticleStatusByIndex(int): Particle index is negative or too big !!" << endl
-         << " It must be inside this range: [0, " << fNParticles-1 << "]" << endl
-         << " Returning null pointer!!" << endl;
+    edm::LogWarning("DatabasePDG")<< "Particle index is negative or too big !!" << endl
+         << " It must be inside this range: (0, " << fNParticles-1 << ")" << endl
+         << " Returning null pointer!!";
     return kFALSE;
   }
   return fStatus[index];
@@ -234,14 +230,14 @@ ParticlePDG* DatabasePDG::GetPDGParticle(int pdg) {
   }
   if(nFindings == 1) return fParticles[firstTimeIndex];
   if(nFindings == 0) {
-    cout << "Warning in DatabasePDG::GetPDGParticle(int): The particle required with PDG = " << pdg
-         << " was not found in the database!!" << endl;
+    edm::LogWarning("DatabasePDG")<< "The particle required with PDG: " << pdg
+         << " was not found in the database!!";
     return 0x0;
   }
   if(nFindings >= 2) {
-    cout << "Warning in DatabasePDG::GetPDGParticle(int): The particle required with PDG = " << pdg
+    edm::LogWarning("DatabasePDG")<< "The particle required with PDG: " << pdg
          << " was found with " << nFindings << " entries in the database. Check it out !!" << endl
-	 << "Returning the first instance found" << endl;
+	 << "Returning the first instance found";
     return fParticles[firstTimeIndex];
   }
   return 0x0;
@@ -258,14 +254,14 @@ bool DatabasePDG::GetPDGParticleStatus(int pdg) {
   }
   if(nFindings == 1) return fStatus[firstTimeIndex];
   if(nFindings == 0) {
-    cout << "Warning in DatabasePDG::GetPDGParticleStatus(int): The particle required with PDG = " << pdg
-         << " was not found in the database!!" << endl;
+    edm::LogWarning("DatabasePDG")<< "The particle required with PDG: " << pdg
+         << " was not found in the database!!";
     return kFALSE;
   }
   if(nFindings >= 2) {
-    cout << "Warning in DatabasePDG::GetPDGParticleStatus(int): The particle status required for PDG = " << pdg
+    edm::LogWarning("DatabasePDG")<< "The particle status required for PDG: " << pdg
          << " was found with " << nFindings << " entries in the database. Check it out !!" << endl
-	 << "Returning the status of first instance found" << endl;
+	 << "Returning the status of first instance found";
     return fStatus[firstTimeIndex];
   }
   return kFALSE;
@@ -282,14 +278,14 @@ ParticlePDG* DatabasePDG::GetPDGParticle(char* name) {
   }
   if(nFindings == 1) return fParticles[firstTimeIndex];
   if(nFindings == 0) {
-    cout << "Warning in DatabasePDG::GetPDGParticle(char*): The particle required with name \"" << name
-         << "\" was not found in the database!!" << endl;
+    edm::LogWarning("DatabasePDG")<< "The particle required with name (" << name
+         << ") was not found in the database!!";
     return 0x0;
   }
   if(nFindings >= 2) {
-    cout << "Warning in DatabasePDG::GetPDGParticle(char*): The particle required with name \"" << name
-         << "\" was found with " << nFindings << " entries in the database. Check it out !!" << endl
-	 << "Returning the first instance found" << endl;
+    edm::LogWarning("DatabasePDG")<< "The particle required with name (" << name
+         << ") was found with " << nFindings << " entries in the database. Check it out !!" << endl
+	 << "Returning the first instance found";
     return fParticles[firstTimeIndex];
   }
   return 0x0;
@@ -306,14 +302,14 @@ bool DatabasePDG::GetPDGParticleStatus(char* name) {
   }
   if(nFindings == 1) return fStatus[firstTimeIndex];
   if(nFindings == 0) {
-    cout << "Warning in DatabasePDG::GetPDGParticleStatus(char*): The particle required with name \"" << name
-         << "\" was not found in the database!!" << endl;
+    edm::LogWarning("DatabasePDG")<< "The particle required with name (" << name
+         << ") was not found in the database!!";
     return kFALSE;
   }
   if(nFindings >= 2) {
-    cout << "Warning in DatabasePDG::GetPDGParticleStatus(char*): The particle status required for name \"" << name
-         << "\" was found with " << nFindings << " entries in the database. Check it out !!" << endl
-	 << "Returning the first instance found" << endl;
+    edm::LogWarning("DatabasePDG")<< "The particle status required for name (" << name
+         << ") was found with " << nFindings << " entries in the database. Check it out !!" << endl
+	 << "Returning the first instance found";
     return fStatus[firstTimeIndex];
   }
   return kFALSE;
@@ -542,7 +538,7 @@ void DatabasePDG::SortParticles() {
 
 
   if(fNParticles<2) {
-    cout << "Warning in DatabasePDG::SortParticles() : No particles to sort. Load data first!!" << endl;
+    edm::LogWarning("DatabasePDG")<< "No particles to sort. Load data first!!";
     return;
   }
 
@@ -587,15 +583,15 @@ int DatabasePDG::GetNParticles(bool all) {
 
 void DatabasePDG::UseThisListOfParticles(char *filename, bool exclusive) {
   if(fNParticles<1) {
-    cout << "Error in DatabasePDG::UseThisListOfParticles(char*, bool) : You must load the data before calling this function!!" << endl;
+    edm::LogError("DatabasePDG")<< "You must load the data before calling this function!!";
     return;
   }
 
   ifstream listFile;
   listFile.open(filename);
   if(!listFile) {
-    cout << "ERROR in DatabasePDG::UseThisListOfParticles(char*, bool) : The ASCII file containing the PDG codes list (\""
-         << filename << "\") was not found !! Aborting..." << endl;
+    edm::LogError("DatabasePDG")<< "The ASCII file containing the PDG codes list ("
+         << filename << ") was not found !!";
     return;
   }
 
@@ -609,7 +605,7 @@ void DatabasePDG::UseThisListOfParticles(char *filename, bool exclusive) {
       listFile >> pdg;
     }
     catch (ios::failure const &problem) {
-      cout << problem.what() << endl;
+      LogDebug("DatabasePDG")<< "ios:failure in list file"<<  problem.what();
       break;
     }
     int found = 0;
@@ -620,12 +616,12 @@ void DatabasePDG::UseThisListOfParticles(char *filename, bool exclusive) {
       }
     }
     if(!found) {
-      cout << "Warning in DatabasePDG::UseThisListOfParticles(char*, bool) : The particle with PDG code "
-	   << pdg << " was asked but not found in the database!!" << endl;
+      edm::LogWarning("DatabasePDG")<< "The particle with PDG code "
+	   << pdg << " was asked but not found in the database!!";
     }
     if(found>1) {
-      cout << "Warning in DatabasePDG::UseThisListOfParticles(char*, bool) : The particle with PDG code "
-	   << pdg << " was found more than once in the database!!" << endl;
+      edm::LogWarning("DatabasePDG")<< "The particle with PDG code "
+	   << pdg << " was found more than once in the database!!";
     }
   }
 

GeneratorInterface/Hydjet2Interface/src/DecayChannel.cc

@@ -39,8 +39,7 @@ DecayChannel::DecayChannel(int mother, double branching, int nDaughters, int *da
   fNDaughters = 0;
   for(int i=0; i<nDaughters; i++) {
     if(i >= kMaxDaughters) {
-      cout << "ERROR in DecayChannel explicit constructor: " << endl;
-      cout << "Number of daughters bigger than the maximum allowed one (" << kMaxDaughters << ") !!" << endl;
+      edm::LogError("DecayChannel")<<"From explicit constructor: Number of daughters bigger than the maximum allowed one (" << kMaxDaughters << ") !!";
     }
     fDaughtersPDG[fNDaughters++] = *(daughters+i);
   }
@@ -49,25 +48,22 @@ DecayChannel::DecayChannel(int mother, double branching, int nDaughters, int *da
 void DecayChannel::SetDaughters(int *daughters, int n) {
   for(int i=0; i<n; i++) {
     if(i >= kMaxDaughters) {
-      cout << "ERROR in DecayChannel::SetDaughters() :" << endl;
-      cout << "Number of daughters bigger than the maximum allowed one (" << kMaxDaughters << ") !!" << endl;
+      edm::LogError("DecayChannel")<<"From SetDaughters(): Number of daughters bigger than the maximum allowed one (" << kMaxDaughters << ") !!";
     }
     fDaughtersPDG[fNDaughters++] = *(daughters+i);
   }
 }
 
 void DecayChannel::AddDaughter(int pdg) {
   if(fNDaughters >= kMaxDaughters) {
-    cout << "ERROR in DecayChannel::AddDaughter() :" << endl;
-    cout << "Number of daughters is already >= than the maximum allowed one (" << kMaxDaughters << ") !!" << endl;
+    edm::LogError("DecayChannel")<<"From AddDaughter(): Number of daughters is already >= than the maximum allowed one (" << kMaxDaughters << ") !!";
   }
   fDaughtersPDG[fNDaughters++] = pdg;
 }
 
 int DecayChannel::GetDaughterPDG(int i) {
   if((i >= fNDaughters) || (i<0)) {
-    cout << "ERROR in DecayChannel::GetDaughterPDG() :" << endl;
-    cout << "Daughter index required is too big or less than zero!! There are only " << fNDaughters << " secondaries in this channel !!" << endl;
+    edm::LogError("DecayChannel")<<"From GetDaughterPDG(): Daughter index required is too big or less than zero!! There are only " << fNDaughters << " secondaries in this channel !!";
     return kNonsensePDG;
   }
   return fDaughtersPDG[i];

GeneratorInterface/Hydjet2Interface/src/GrandCanonical.cc

@@ -6,7 +6,6 @@ November. 2, 2005
 
 */
 #include <iostream>
-#include <TError.h>
 #include <TMath.h>
 
 #include "GeneratorInterface/Hydjet2Interface/interface/GrandCanonical.h"