---

CellMaker.cxx

Go to the documentation of this file.

// ================================================
// CellMaker 
// ================================================
//
// THIS TEXT TO BE REPLACED BY ATLAS STANDARD FORMAT
//
// Namespace Atlfast
//

#include "AtlfastAlgs/CellMaker.h"
#include "AtlfastAlgs/CellsAboveThreshold.h"
#include "AtlfastAlgs/TransportedParticle.h"
#include "AtlfastAlgs/TransportedParticleCollection.h"
#include "AtlfastAlgs/Calorimeter.h"
#include "AtlfastAlgs/MagField.h"
#include "AtlfastAlgs/ICellSelector.h"
#include "AtlfastAlgs/ISmearer.h"
#include "AtlfastAlgs/CellSmearer.h"
#include "AtlfastAlgs/GlobalEventData.h"
#include "AtlfastAlgs/EPileupMap.h"

#include "AtlfastEvent/CollectionDefs.h"

#include "AtlfastUtils/TesIO.h"
#include "AtlfastUtils/HeaderPrinter.h"

#include <cmath>
#include <iomanip>
#include <vector>
// Gaudi includes
#include "GaudiKernel/DataSvc.h"
#include "GaudiKernel/Chrono.h"
#include "StoreGate/DataHandle.h"
#include "GaudiKernel/ISvcLocator.h"
#include "GaudiKernel/MsgStream.h"
//Generator includes
//#include "GeneratorObjects/McEventCollection.h"

namespace Atlfast {
  CellMaker:: CellMaker(const std::string& name, ISvcLocator* pSvcLocator):
    Algorithm(name,pSvcLocator){
    
    //default paths for entities in the TES
    m_inputLocation       = "/Event/McEventCollection";
    m_outputLocation      = "/Event/AtlfastCells" ;
    
    // default parameters to build the calorimeter
    m_etaCoverage        = 5.0;
    m_minETCell          = 0.0;
    m_granBarrelEta      = 0.1;
    m_granBarrelPhi      = 0.1;
    m_granForwardEta     = 0.2;
    m_granForwardPhi     = 0.2;
    m_doSmearing         = false;
    m_fastShower         = false;
    
    declareProperty("InputLocation",      m_inputLocation);
    declareProperty("OutputLocation",     m_outputLocation);
    declareProperty("EtaCoverage",        m_etaCoverage);
    declareProperty("MinETCell",          m_minETCell);
    declareProperty("GranBarrelEta",      m_granBarrelEta);
    declareProperty("GranBarrelPhi",      m_granBarrelPhi);
    declareProperty("GranForwardEta",     m_granForwardEta);
    declareProperty("GranForwardPhi",     m_granForwardPhi);
    declareProperty("DoSmearing",     m_doSmearing);
    declareProperty("FastShower",     m_fastShower);
    
  }
  
  //__________________________________________________________________________
  CellMaker::~CellMaker(){

    MsgStream log( messageService(), name() );
    log<<MSG::DEBUG<<"CellMaker destructor starts "<<endreq;

    if (m_calorimeter) {
      delete m_calorimeter;
    }
    if (m_tesIO) {
      delete m_tesIO;
    }
    if (m_cellSelector) {
      delete m_cellSelector;
    }
    if (m_magField) {
      delete m_magField;
    }
    log<<MSG::DEBUG<<"CellMaker destructor ends "<<endreq;
  }
  //__________________________________________________________________________
  StatusCode CellMaker::initialize()
  {
    
    MsgStream log( messageService(), name() );

    //set up objects used to determine particle charge,
    // and to select input MC 4-vectors.
    m_cellSelector  = new CellsAboveThreshold(m_minETCell);
    
    //    m_tesIO = new TesIO(eventDataService());
    m_tesIO = new TesIO();

    //get the Global Event Data using singleton pattern
    GlobalEventData* ged = GlobalEventData::Instance();
    m_fieldOn          = ged->fieldOn();
    m_mcSelector       = ged -> visibleToCal();
    m_barrelForwardEta = ged -> barrelForwardEta();
    m_magField      = new MagField(m_fieldOn);
    if (fabs(m_etaCoverage) < fabs(m_barrelForwardEta) ) {
      log << MSG::ERROR 
          << "asked to makecalo with max extent " 
          << m_etaCoverage 
          << "but barrel-forward boundary at " 
          << m_barrelForwardEta << endreq;
      return StatusCode::FAILURE;
    }
    m_calorimeter = new Calorimeter(
                                    log,
                                    m_fastShower,
                                    m_etaCoverage, 
                                    m_barrelForwardEta,
                                    m_granBarrelEta, 
                                    m_granBarrelPhi, 
                                    m_granForwardEta, 
                                    m_granForwardPhi
                                    );
    if(m_doSmearing == true){
      m_epileupMap = new EPileupMap(ged->lumi());
      m_smearer = new CellSmearer(ged->randSeed(),m_barrelForwardEta);
    }else{
      m_epileupMap = NULL;
      m_smearer = NULL;
    }

    HeaderPrinter hp("Atlfast Cell Maker:", log);
    
    hp.add("min Cell ET             ", m_minETCell);  
    hp.add("Magnetic Field is on?   ", m_fieldOn);
    hp.add("Input Location          ", m_inputLocation);    
    hp.add("Output Location         ", m_outputLocation);
    hp.add(" DoSmearing             ", m_doSmearing);    
    hp.add(" Run FastShower         ", m_fastShower);    
    hp.add("Calorimeter Geometry:   ");
    hp.add(" Total Eta Range        ", m_etaCoverage);
    hp.add(" Barrel Eta Range       ", m_barrelForwardEta);
    hp.add(" Barrel Eta Granularity ", m_granBarrelEta);
    hp.add(" Barrel Phi Granularity ", m_granBarrelPhi);
    hp.add(" Forward Eta Granularity", m_granForwardEta);
    hp.add(" Forward Phi Granularity", m_granForwardPhi);
    hp.print();


      
    return StatusCode::SUCCESS; 
  }
  
  StatusCode CellMaker::finalize()
  {
    
    MsgStream log( messageService(), name() );
    log<<MSG::INFO<<"Finalizing"<<endreq;
    return StatusCode::SUCCESS; 
  }
  
  //_________________________________________________________________________
  StatusCode CellMaker::execute()
  {
    //.............................................
    
    MsgStream log( messageService(), name() );
    log << MSG::DEBUG << "CellMaker execute()" << endreq;

    std::string mess;
    

    // read MC particles from TES
    MCparticleCollection  mcParticles ;
    TesIoStat stat = m_tesIO->getMC( mcParticles, m_mcSelector ) ;
    mess = stat? "Retrieved MC from TES ":"Failed MC retrieve from TES";
    log << MSG::DEBUG << mess << endreq;
    
    // calculate phi after the bending in the magnetic field, 
    TransportedParticleCollection particlesAtCal;
    (*m_magField)(mcParticles, particlesAtCal);
    log << MSG::DEBUG << "bend" << endreq;
    
    //    return stat;
    // give to calorimeter (iseq gives the interval [begin,end)
    TransportedParticleCollectionIter fpart = particlesAtCal.begin();
    TransportedParticleCollectionIter lpart = particlesAtCal.end();
    m_calorimeter->deposit(fpart, lpart);


    // Add in the energy pileup 
    if(m_doSmearing){
      //Smear the cells resolution first
      m_calorimeter->smearCells(m_smearer);
      log << MSG::DEBUG << "Doing Smearing" << endreq;
      m_epileupMap->fillMap();
      log << MSG::DEBUG 
          << "Size of EM EPileup vector: " 
          << m_epileupMap->getEMMap()->size() 
          << endreq;
      std::vector<Atlfast::EPileupDeposit*>::iterator fdep = 
        m_epileupMap->getEMMap()->begin();
      std::vector<Atlfast::EPileupDeposit*>::iterator ldep = 
        m_epileupMap->getEMMap()->end();
      m_calorimeter->deposit(fdep, ldep);
      log << MSG::DEBUG << "ECAL pileup deposited" << endreq;
      log << MSG::DEBUG << "Size of HAD EPileup vector: " 
          << m_epileupMap->getHadMap()->size() 
          << endreq;
      fdep = m_epileupMap->getHadMap()->begin();
      ldep = m_epileupMap->getHadMap()->end();
      m_calorimeter->deposit(fdep, ldep);
      log << MSG::DEBUG << "HCAL pileup deposited" << endreq;
    }

    log << MSG::DEBUG << "deposit" << endreq;
    
    // ask for Cells passing selection criteria of m_cellSelector (=Et).
    //************************
    //this worked 23/12 PS
    ITwoCptCellCollection* cells = new ITwoCptCellCollection;
    //************************
    m_calorimeter->giveHitCells(m_cellSelector, cells);
    log << MSG::DEBUG << "hit cells retrieved, " << cells->size() 
        << " cells were hit this event " << endreq;

    //    return stat;

    // finished with the calorimeter
    m_calorimeter->reset();
    log << MSG::DEBUG << "calorimeter reset" << endreq;
    //    return stat;
    //    return stat;
    
    // then store and return
    stat = m_tesIO->store(cells, m_outputLocation );  
    mess = stat? " Stored cells to TES":" Failed store cells to TES";
    log << MSG::DEBUG << mess << endreq;

    return stat;
  }

} // end of namespace bracket

---