---

digiRootReaderAlg.cxx

Go to the documentation of this file.

#include "GaudiKernel/MsgStream.h"
#include "GaudiKernel/AlgFactory.h"
#include "GaudiKernel/IDataProviderSvc.h"
#include "GaudiKernel/SmartDataPtr.h"
#include "GaudiKernel/Property.h"
#include "GlastSvc/GlastDetSvc/IGlastDetSvc.h"

#include "GaudiKernel/DataObject.h"
#include "GaudiKernel/DataSvc.h"
#include "GaudiKernel/ObjectVector.h"


#include "digiRootReaderAlg.h"
// switch to new root Digi package
#include "digiRootData/AcdTile.h"
#include "digiRootData/CalLog.h"
#include "digiRootData/StripId.h"
#include "digiRootData/TkrLayer.h"
//#include "digiRootData/Tagger.h"

#include "TROOT.h"
#include "TFile.h"
#include "TTree.h"
#include "TDirectory.h"

#include "digiRootData/DigiEvent.h"

#include "GlastEvent/Hits/SiLayers.h"
#include "GlastEvent/TopLevel/Event.h"
#include "CalRecon/CalADCLogs.h"

//#include "GlastEvent/Digi/AcdDigi.h"

static const AlgFactory<digiRootReaderAlg>  Factory;
const IAlgFactory& digiRootReaderAlgFactory = Factory;


digiRootReaderAlg::digiRootReaderAlg(const std::string& name, ISvcLocator* pSvcLocator) :
Algorithm(name, pSvcLocator) 
{
    declareProperty("rawRootFile",m_fileName);
    m_numEvents = 0;
    m_currentEvent = 0;
    m_eventId = 0;
    m_runId = 0;
}


StatusCode digiRootReaderAlg::initialize()
{
    StatusCode sc = StatusCode::SUCCESS;
    MsgStream log(msgSvc(), name());
    // Use the Job options service to set the Algorithm's parameters
    setProperties();
    
    IProperty* glastPropMgr=0;
    sc = serviceLocator()->getService("EventSelector", IID_IProperty,
                         reinterpret_cast<IInterface*&>( glastPropMgr ));
    if( sc.isFailure() ) return sc;
      
    IntegerProperty evtMax("EvtMax",0);
    sc = glastPropMgr->getProperty( &evtMax );
    if (sc.isFailure()) return sc;

    m_maxEvent = evtMax.value();


    if(0 == gROOT)
    {
        static TROOT root("root","Root I/O");
        log << MSG::INFO << "Initializing ROOT" << endreq;
    } else {
        log << MSG::INFO << "Root already initialized(not an error)" <<endreq;
    }
    
    // Setup the input raw root file for reading
    TDirectory *saveDir = gDirectory;
    m_rawFile = new TFile(m_fileName.c_str(), "READ");
    if (!(m_rawFile->IsOpen())) return StatusCode::FAILURE;

    m_rawTree = (TTree*) m_rawFile->Get("T");
    if (!m_rawTree) return StatusCode::FAILURE;
    
    m_numEvents = m_rawTree->GetEntries();
    
    m_event = new DigiEvent();
    m_event->Clean();
    
    m_rawTree->SetBranchAddress("DigiEvent", &m_event);
    
    
    saveDir->cd();
    return sc;
}

StatusCode digiRootReaderAlg::execute()
{
    StatusCode sc = StatusCode::SUCCESS;
    MsgStream   log( msgSvc(), name() );
    
    if (m_currentEvent >= m_numEvents) {
        log << MSG::INFO << "reached end of raw root file " << endreq;
        return StatusCode::FAILURE;
    }
    
    if (m_currentEvent >= m_maxEvent) {
        log << MSG::INFO << "Finished processing the requested number of evnets " << endreq;
        return StatusCode::FAILURE;
    }

    readEvent();
   
    sc = retrieve();
    load(m_SiLayers);
    load(m_CalADCLogs);
    load(m_AcdDigiList);
    
    sc = storeHeader();

    ++m_currentEvent;
    
    return sc;
}

void digiRootReaderAlg::readEvent() {
    
    TDirectory *saveDir = gDirectory;
    m_event->Clean();   
    m_rawFile->cd();
    
    Int_t nb = 0,nbytes = 0;
    nb = m_rawTree->GetEvent((Int_t) m_currentEvent);
    nbytes += nb;
    
    m_eventId = m_event->getEventId();
    m_runId = m_event->getRunId();

    saveDir->cd();
    
    return;
}

StatusCode digiRootReaderAlg::storeHeader() {
    MsgStream log(msgSvc(), name());
    StatusCode sc; 

    // There is no proper Event object in the TDS - since
    // we are not using the EventSelector.  Instead a generic
    // Data Object is the root of the TDS.  In this case, we
    // create our own Event object and store it in the TDS.
    Event *eventHeader = new Event();
 
    eventHeader->setEvent(m_eventId);
    eventHeader->setRun(m_runId);
    
    sc = eventSvc()->registerObject("/Event/Header", eventHeader);        
    
    return sc;
}

StatusCode digiRootReaderAlg::finalize()
{
    close();
    
    StatusCode sc = StatusCode::SUCCESS;
    return sc;
}

void digiRootReaderAlg::close() 
{
    TDirectory *saveDir = gDirectory;
    m_rawFile->cd();
    m_rawFile->Close();
    saveDir->cd();
}

StatusCode digiRootReaderAlg::retrieve()
{
    MsgStream log(msgSvc(), name());
    StatusCode sc = StatusCode::SUCCESS;
    
    m_SiLayers = 0;
    m_SiLayers = new SiLayers();
    m_SiLayers->clear();
    
    
    m_CalADCLogs = 0;
    m_CalADCLogs = new CalADCLogs();

    m_AcdDigiList = new AcdDigiVector;
    
    
    DataObject* pnode=0;
    
    sc = eventSvc()->retrieveObject( "/Event/TkrRecon", pnode );
    
    if( sc.isFailure() ) {
        sc = eventSvc()->registerObject("/Event/TkrRecon",new DataObject);
        if( sc.isFailure() ) {
            
            log << MSG::ERROR << "Could not create TkrRecon directory" << endreq;
            return sc;
        }
    }
    
    
    StatusCode sc_tkr = eventSvc()->registerObject("/Event/TkrRecon/SiLayers",m_SiLayers);
    
    sc = eventSvc()->retrieveObject( "/Event/CalRecon", pnode );
    
    if( sc.isFailure() ) {
        sc = eventSvc()->registerObject("/Event/CalRecon",new DataObject);
        if( sc.isFailure() ) {
            
            log << MSG::ERROR << "Could not create CalRecon directory" << endreq;
            return sc;
        }
    }
    
    StatusCode sc_cal = eventSvc()->registerObject("/Event/CalRecon/CalADCLogs",m_CalADCLogs);    
    
    sc = eventSvc()->retrieveObject("/Event/Digi", pnode);
    if (sc.isFailure()) {
        sc = eventSvc()->registerObject("/Event/Digi", new DataObject);
        if (sc.isFailure()) {
            log << MSG::ERROR << "Could not create Digi directory" << endreq;
            return sc;
        }
    }
    
    
    sc = eventSvc()->registerObject("/Event/Digi/AcdDigis", m_AcdDigiList);
    if (sc.isFailure()) {
        log << MSG::ERROR << "could not load AcdDigs" << endreq;
    }

    if (sc_cal.isFailure() || sc_tkr.isFailure()) sc = StatusCode::FAILURE;
    
    return sc;
}

void digiRootReaderAlg::load(SiLayers* siLayers)
{
    int nTKR = m_event->getTkrDigi()->GetEntries();
    for (int i = 0; i < nTKR; i++) {
        
        TkrLayer* rLayer = (TkrLayer*) m_event->getTkrDigi()->At(i);
        // For the testbeam ROOT classes.. TkrLayer::X was defined 
        // as a layer ORIENTED in the X direction
        // This has been fixed... TkrLayer::X now means MEASURES X
        // So we replace the following line:
        // int XY = (rLayer->getXY() == TkrLayer::X ? 1 : 0);
        // with this new line:
        int XY = rLayer->getXY();

        double tot = rLayer->getToT(0);
        
        int nStrips = rLayer->getStrips()->GetEntries();
        if (nStrips >0 ) {
            SiLayer* siLayer = new SiLayer(rLayer->getPlaneNum(), XY, tot);
            StripId* rStrip;
            for (int iStrip=0; iStrip < nStrips; iStrip++) {
                rStrip = (StripId*) rLayer->getStrips()->At(iStrip);
                int iaddress = rStrip->getId(); 
                siLayer->addStrip((int) iaddress);
            }
            siLayers->add(siLayer);
        }
    }
}



void digiRootReaderAlg::load(CalADCLogs* calLogs)
{       
    int iread = 0;
    double BASE = 4095; // the 0 for the ADC's
    int NLAYERS = 3;
    const MAXLAYERNUM = 7;

    int nCAL = m_event->getCalDigi()->GetEntries();
    for (int i = 0; i < nCAL; i++) {
        CalLog* rLog = (CalLog*) m_event->getCalDigi()->At(i);        
        LogId* rLogID= rLog->getLogId();
      
        // The layer numbers in the Raw Root file adhere to the
        // Ritz convention - layer num == 0 is the front most layer
        // closest to the TKR.  Here, iplane is numbered from the bottom (0)
        // and is [0,3] and where iview is used to denote the orientation 
        // of the layer
        int iplane = (MAXLAYERNUM - rLogID->getLayer()) / 2;
        // HMK:  The columns ids in the Raw Root file now conform to the standard numbering scheme
        // I believe, that CalRecon also adheres to the standard in the respect
        // So there is no need to modify the column number obtained from the raw root file
        int icolumn = rLogID->getColumn();
        // Here, the Raw Root file provides getXY == 0 to mean a log oriented along the 
        // X-axis, so it measures Y.  HMK believes that iview should be 1 in this case
        // and 0 for logs oriented along the Y-axis.
        int iview = 1-(rLogID->getXY());

        // Old Code that will no longer work
        /*
        int ADCID   = rLogID->getId();
        int iplane = ADCID/2;
        iplane = NLAYERS - (iplane%4);
        int icolumn = ADCID/8;
        int iview   = 1-ADCID%2;  
        */

        // We have fixed the digi Root classes so that the 
        // columns now conform to the Steve Ritz numbering standard
        // No need to change the orientation.
        //if (iview == 0) icolumn = 9 - icolumn; 
        
        CalADCLog* log = calLogs->getLogID(CalLogID::ID(iplane,CalDetGeo::makeAxis(iview),icolumn));
        
        log->setNegADC(CalBase::LEX,BASE-rLog->getAdc(CalLog::NEG,CalLog::LowEnX4));
        log->setNegADC(CalBase::LE, BASE-rLog->getAdc(CalLog::NEG,CalLog::LowEnX1));
        log->setNegADC(CalBase::HEX,BASE-rLog->getAdc(CalLog::NEG,CalLog::HiEnX8));
        log->setNegADC(CalBase::HE, BASE-rLog->getAdc(CalLog::NEG,CalLog::HiEnX1));
        
        log->setPosADC(CalBase::LEX,BASE-rLog->getAdc(CalLog::POS,CalLog::LowEnX4));
        log->setPosADC(CalBase::LE, BASE-rLog->getAdc(CalLog::POS,CalLog::LowEnX1));
        log->setPosADC(CalBase::HEX,BASE-rLog->getAdc(CalLog::POS,CalLog::HiEnX8));
        log->setPosADC(CalBase::HE, BASE-rLog->getAdc(CalLog::POS,CalLog::HiEnX1));
    }
}


void digiRootReaderAlg::load(AcdDigiVector *acdDigiVec) {

    int numACDTiles = m_event->getAcdDigi()->GetEntries();
    for (int i = 0; i < numACDTiles; i++) {
        AcdTile* tile = (AcdTile*) m_event->getAcdDigi()->At(i);        
        AcdId* id= tile->getId();

        acdDigiVec->push_back(new AcdDigi(idents::AcdId(id->getId(2)), tile->getPulseHeight(), 
            tile->getVeto(), false, false) );
    }

}

---