add HI printouts (TruthIO-00-00-13)

        * add HI printouts in PrintMC (A. Olszewski)
        * tag as TruthIO-00-00-13

2015-05-05 E. M. Lobodzinska <ewelina@mail.desy.de>
        * inherit ReadHepEvtFromAscii from AthAlgorithm
        * tag as TruthIO-00-00-12

2015-04-17 E. M. Lobodzinska <ewelina@mail.desy.de>
        * add ReadHepEvtFromAscii (Dirk Sammel)
        * tag as TruthIO-00-00-11

Generators/TruthIO/TruthIO/ReadHepEvtFromAscii.h

+/*
+  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "AthenaBaseComps/AthAlgorithm.h"
+#include "HepMC/IO_GenEvent.h"
+
+class StoreGateSvc;
+
+class ReadHepEvtFromAscii:public AthAlgorithm {
+public:
+  ReadHepEvtFromAscii(const std::string& name, ISvcLocator* pSvcLocator);
+  StatusCode initialize();
+  StatusCode execute();
+  StatusCode finalize();
+
+private:
+
+  StoreGateSvc* m_sgSvc;
+  HepMC::IO_GenEvent* m_ascii_in;
+  
+  // Setable Properties:-
+  std::string m_key; 
+  std::string m_input_file;
+  
+};

Generators/TruthIO/cmt/requirements

@@ -6,8 +6,10 @@ use AtlasPolicy      AtlasPolicy-*
 use GaudiInterface   GaudiInterface-*      External
 use GeneratorModules GeneratorModules-*    Generators
 use AtlasHepMC	     AtlasHepMC-*	       External
+use AthenaBaseComps     AthenaBaseComps-*      Control
 
 private
+use StoreGate        StoreGate-*           Control
 use GeneratorObjects GeneratorObjects-*    Generators
 use EventInfo        EventInfo-*           Event
 use HepPDT           v*                    LCG_Interfaces
@@ -15,10 +17,11 @@ end_private
 
 private
 #macro_append HepMC_linkopts " $(HepMC_IO_linkopts) "
+macro_append  TruthIO_shlibflags  "$(HepMC_IO_linkopts)"
 end_private
 
 public
-library TruthIO *.cxx -s=components *.cxx
+library TruthIO *.cxx *.F -s=components *.cxx
 apply_pattern component_library
 
 apply_pattern declare_joboptions files="*.py"

Generators/TruthIO/src/PrintMC.cxx

@@ -105,6 +105,27 @@ StatusCode PrintMC::execute() {
       std::cout << " Entries this event: " << evt->vertices_size() << " vertices, "
                 << evt->particles_size() << " particles.\n";
 
+      if (evt->heavy_ion()) {
+        std::cout << " HeavyIon: jatt=" << evt->heavy_ion()->Ncoll_hard()
+                  << " np=" << evt->heavy_ion()->Npart_proj()
+                  << " nt=" << evt->heavy_ion()->Npart_targ()
+                  << " ncoll=" << evt->heavy_ion()->Ncoll()
+                  << " specn=" << evt->heavy_ion()->spectator_neutrons()
+                  << " specp=" << evt->heavy_ion()->spectator_protons()
+                  << " n01=" << evt->heavy_ion()->N_Nwounded_collisions()
+                  << " n10=" << evt->heavy_ion()->Nwounded_N_collisions()
+                  << " n11=" << evt->heavy_ion()->Nwounded_Nwounded_collisions()
+                  << " impact=" << evt->heavy_ion()->impact_parameter()
+                  << " evplane=" << evt->heavy_ion()->event_plane_angle()
+                  << " ecc=" << evt->heavy_ion()->eccentricity()
+                  << " sigmaNNinel=" << evt->heavy_ion()->sigma_inel_NN()
+		  << std::endl;
+      }
+      else {
+        std::cout << "HeavyIon: EMPTY"
+		  << std::endl;
+      }
+
       // Random State
       //	std::cout << " RndmState(" << evt->random_states().size() << ")=";
       //for ( std::vector<long int>::const_iterator rs
@@ -132,7 +153,8 @@ StatusCode PrintMC::execute() {
 		  << std::endl;
       }
       else {
-        std::cout << "PdfInfo: EMPTY";
+        std::cout << "PdfInfo: EMPTY"
+		  << std::endl;
       }
 
       // Print a legend to describe the particle info

Generators/TruthIO/src/ReadHepEvtFromAscii.cxx

+/*
+  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include <set>
+
+#include "TruthIO/ReadHepEvtFromAscii.h"
+#include "GeneratorObjects/McEventCollection.h"
+#include "HepMC/GenEvent.h"
+#include "HepMC/IO_HEPEVT.h"
+#include "HepMC/HEPEVT_Wrapper.h"
+
+#include "GaudiKernel/DataSvc.h"
+
+#include "StoreGate/StoreGateSvc.h"
+
+extern "C" {
+  void openfile_(int & ifile, const char* filename, int);
+  void readevt_(int& ifile);
+  void closefile_(int& ifile);
+}
+
+
+ReadHepEvtFromAscii::ReadHepEvtFromAscii(const std::string& name, ISvcLocator* pSvcLocator) 
+  :  AthAlgorithm(name, pSvcLocator),
+     m_sgSvc(0), m_ascii_in(0)
+{
+  // Set users' request
+  declareProperty("McEventKey",  m_key        = "GEN_EVENT");
+  declareProperty("HepEvtFile",   m_input_file = "events.ascii");
+}
+
+StatusCode ReadHepEvtFromAscii::initialize(){
+
+  msg(MSG::INFO) << ">>> ReadHepEvtFromAscii from Initialize" << endreq;
+
+  StatusCode sc = service("StoreGateSvc", m_sgSvc);
+    
+  if (sc.isFailure()) {
+    msg(MSG::ERROR) << "Could not find StoreGateSvc" << endreq;
+    return sc;
+  }
+
+  HepMC::HEPEVT_Wrapper::set_sizeof_int(4);
+  HepMC::HEPEVT_Wrapper::set_sizeof_real(8);
+  HepMC::HEPEVT_Wrapper::set_max_number_entries(10000);
+  
+  // Initialize input file
+//  m_ascii_in = new HepMC::IO_Ascii(m_input_file.c_str(), std::ios::in);
+  int ifile=5;
+  closefile_(ifile);
+  openfile_(ifile,m_input_file.c_str(),m_input_file.size());
+
+  // Initialization terminated
+  return StatusCode::SUCCESS;
+}
+StatusCode ReadHepEvtFromAscii::execute() {
+
+  msg(MSG::INFO) << ">>> ReadHepEvtFromAscii from execute" << endreq;
+  
+  McEventCollection* mcEvtColl;
+
+  if ( m_sgSvc->contains<McEventCollection>(m_key) && m_sgSvc->retrieve(mcEvtColl, m_key).isSuccess() ) {
+    if (msgLvl(MSG::VERBOSE)) msg(MSG::VERBOSE) << "found an McEventCollecion in store" << endreq;
+  } else {
+    // McCollection doesn't exist. Create it (empty)
+    if (msgLvl(MSG::VERBOSE)) msg(MSG::VERBOSE) << "create new McEventCollecion in store" << endreq;
+    mcEvtColl = new McEventCollection;
+    StatusCode status = m_sgSvc->record( mcEvtColl, m_key );
+    if (status.isFailure()) {
+      msg(MSG::ERROR) << "Could not record McEventCollection" << endreq;
+      return status;
+    }
+  }
+
+  HepMC::GenEvent* evt = new HepMC::GenEvent(); 
+  int ifile=5; 
+  readevt_(ifile);
+
+  HepMC::IO_HEPEVT hepio;
+  hepio.set_print_inconsistency_errors(0);
+  hepio.fill_next_event(evt);
+ 
+  // Convert GeV to MeV 
+  for ( HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p )
+      {
+	  HepMC::FourVector newMomentum(0.,0.,0.,0.);
+	  newMomentum.setPx( (*p)->momentum().px() * 1000. );
+	  newMomentum.setPy( (*p)->momentum().py() * 1000. );
+	  newMomentum.setPz( (*p)->momentum().pz() * 1000. );
+	  newMomentum.setE( (*p)->momentum().e() * 1000. );
+	  (*p)->set_momentum(newMomentum);
+
+      }
+  mcEvtColl->push_back(evt);
+  
+  // End of execution for each event
+  return StatusCode::SUCCESS;
+}
+StatusCode ReadHepEvtFromAscii::finalize() {
+
+  msg(MSG::INFO) << ">>> ReadHepEvtFromAscii from finalize" << endreq;
+
+  int ifile=5;
+  closefile_(ifile);
+  // End of finalization step
+  return StatusCode::SUCCESS;
+
+}

Generators/TruthIO/src/components/TruthIO_entries.cxx

@@ -4,12 +4,14 @@
 #include "TruthIO/WriteHepMC.h"
 #include "TruthIO/PrintMC.h"
 #include "TruthIO/PrintHijingPars.h"
+#include "TruthIO/ReadHepEvtFromAscii.h"
 
 DECLARE_ALGORITHM_FACTORY(DumpMC)
 DECLARE_ALGORITHM_FACTORY(HepMCReadFromFile)
 DECLARE_ALGORITHM_FACTORY(PrintMC)
 DECLARE_ALGORITHM_FACTORY(WriteHepMC)
 DECLARE_ALGORITHM_FACTORY(PrintHijingPars)
+DECLARE_ALGORITHM_FACTORY(ReadHepEvtFromAscii)
 
 DECLARE_FACTORY_ENTRIES(GeneratorModules) {
   DECLARE_ALGORITHM(DumpMC)
@@ -17,4 +19,5 @@ DECLARE_FACTORY_ENTRIES(GeneratorModules) {
   DECLARE_ALGORITHM(PrintMC)
   DECLARE_ALGORITHM(WriteHepMC)
   DECLARE_ALGORITHM(PrintHijingPars)
+  DECLARE_ALGORITHM(ReadHepEvtFromAscii)
 }

Generators/TruthIO/src/readEvt.F

+      subroutine openFile(ifile,filename)
+      integer ifile
+      character*(*) filename
+      logical ex
+
+      inquire(FILE=filename,EXIST=ex)
+      if(.not. ex) write(*,*) 'File does not exist:',filename
+      open(ifile,FILE=filename,status='old')
+      end
+
+      subroutine readEvt(ifile)
+      integer ifile
+c...Hepevt common block
+      integer nmxhep,nevhep,nhep,isthep,idhep,jmohep,jdahep
+      double precision phep,vhep
+      PARAMETER (NMXHEP=10000)
+      COMMON/HEPEVT/NEVHEP,NHEP,ISTHEP(NMXHEP),IDHEP(NMXHEP),
+     & JMOHEP(2,NMXHEP),JDAHEP(2,NMXHEP),PHEP(5,NMXHEP),VHEP(4,NMXHEP)
+
+      integer i
+
+      read(ifile,*) nevhep,nhep
+      do i=1,nhep
+        read(ifile,*) dummy,isthep(i),idhep(i),jmohep(1,i),
+     & jmohep(2,i),jdahep(1,i),jdahep(2,i)
+        read(ifile,*) phep(1,i),phep(2,i),phep(3,i),phep(4,i),phep(5,i)
+        read(ifile,*) vhep(1,i),vhep(2,i),vhep(3,i),vhep(4,i)
+      enddo
+
+      end
+
+      subroutine closeFile(ifile)
+      integer ifile
+
+      close(ifile)
+
+      end