RecMdcDedxCnv.cxx

Go to the documentation of this file.

#ifndef RecMdcDedxCnv_CXX
#define RecMdcDedxCnv_CXX 1
#include "EventModel/EventModel.h"
#include "GaudiKernel/DataObject.h"
#include "GaudiKernel/MsgStream.h"
#include "GaudiKernel/ObjectVector.h"
#include "MdcRecEvent/RecMdcDedx.h"
#include "MdcRecEvent/RecMdcDedxHit.h"
#include "MdcRecEvent/RecMdcKalTrack.h"
#include "MdcRecEvent/RecMdcTrack.h"
#include "ReconEvent/ReconEvent.h"
#include "TClonesArray.h"
 
#include "RootEventData/TRecMdcDedx.h" // standard root object
#include "RootEventData/TRecTrackEvent.h"
 
#include "RootCnvSvc/Rec/RecMdcDedxCnv.h"
#include "RootCnvSvc/Rec/RecTrackCnv.h"
#include "RootCnvSvc/RootAddress.h"
 
#include <vector>
 
using namespace std;
 
// Instantiation of a static factory class used by clients to create
// instances of this service
// static CnvFactory<RecMdcDedxCnv> s_factory;
// const ICnvFactory& RecMdcDedxCnvFactory = s_factory;
 
RecMdcDedxCnv::RecMdcDedxCnv( ISvcLocator* svc ) : RootEventBaseCnv( classID(), svc ) {
 
  // Here we associate this converter with the /Event path on the TDS.
  MsgStream log( msgSvc(), "RecMdcDedxCnv" );
  // log << MSG::DEBUG << "Constructor called for " << objType() << endmsg;
  // m_rootTreename ="Rec";
  m_rootBranchname = "m_recMdcDedxCol";
  // declareObject(EventModel::Recon::RecMdcDedxCol, objType(), m_rootTreename,
  // m_rootBranchname);
  m_adresses.push_back( &m_recMdcDedxCol );
  m_recMdcDedxCol = 0;
}
 
StatusCode RecMdcDedxCnv::TObjectToDataObject( DataObject*& refpObject ) {
  // creation of TDS object from root object
  MsgStream log( msgSvc(), "RecMdcDedxCnv" );
  log << MSG::DEBUG << "RecMdcDedxCnv::TObjectToDataObject" << endmsg;
  StatusCode sc = StatusCode::SUCCESS;
 
  // create the TDS location for the Dedx Collection
  RecMdcDedxCol* recMdcDedxCol = new RecMdcDedxCol;
  refpObject                   = recMdcDedxCol;
 
  if ( !m_recMdcDedxCol ) return sc;
 
  //*******************Dst to Rec convert in case of no Rec data in
  // TDS****************************
  IDataProviderSvc* dataSvc = 0;
  sc = serviceLocator()->getService( "EventDataSvc", IDataProviderSvc::interfaceID(),
                                     (IInterface*&)dataSvc );
  if ( sc.isFailure() )
  {
    log << MSG::FATAL << "Could not get EventDataSvc in RecMdcDedxCnv" << endmsg;
    return sc;
  }
 
  //----------get hit smartRefVector--------
  SmartDataPtr<RecMdcDedxHitCol> recMdcDedxHitCol( dataSvc, "/Event/Recon/RecMdcDedxHitCol" );
  if ( !recMdcDedxHitCol )
  {
    log << MSG::FATAL << "Could not find RecMdcDedxHitCol" << endmsg;
    return ( StatusCode::FAILURE );
  }
 
  int trackID;
  // Retrieve Mdc Track
  SmartDataPtr<RecMdcTrackCol> mdcTrackCol( dataSvc, EventModel::Recon::RecMdcTrackCol );
  if ( !mdcTrackCol )
  {
    log << MSG::INFO << "Could not find RecMdcTrackCol" << endmsg;
    SmartDataPtr<DstMdcTrackCol> dstMdcTrackCol( dataSvc, "/Event/Dst/DstMdcTrackCol" );
    if ( !dstMdcTrackCol ) { log << MSG::INFO << "Could not find DstMdcTrackCol" << endmsg; }
    else
    {
      RecMdcTrackCol*          mdcTrackCol = new RecMdcTrackCol();
      DstMdcTrackCol::iterator iter_mdc    = dstMdcTrackCol->begin();
      trackID                              = 0;
      for ( ; iter_mdc != dstMdcTrackCol->end(); iter_mdc++, trackID++ )
      {
        RecMdcTrack* recMdcTrack = new RecMdcTrack();
        *recMdcTrack             = **iter_mdc;
        ( *mdcTrackCol ).push_back( recMdcTrack );
        log << MSG::INFO << " Mdc Track ID = " << trackID
            << " Mdc Track Nster = " << ( *iter_mdc )->nster() << endmsg;
      }
      sc = dataSvc->registerObject( EventModel::Recon::RecMdcTrackCol, mdcTrackCol );
    }
  }
 
  // Retrieve MdcKal track
  SmartDataPtr<RecMdcKalTrackCol> mdcKalTrackCol( dataSvc,
                                                  EventModel::Recon::RecMdcKalTrackCol );
  if ( !mdcKalTrackCol )
  {
    log << MSG::INFO << "Could not find RecMdcKalTrackCol" << endmsg;
    SmartDataPtr<DstMdcKalTrackCol> dstMdcKalTrackCol( dataSvc,
                                                       "/Event/Dst/DstMdcKalTrackCol" );
    if ( !dstMdcKalTrackCol )
    { log << MSG::INFO << "Could not find DstMdcKalTrackCol" << endmsg; }
    else
    {
      RecMdcKalTrackCol*          mdcKalTrackCol = new RecMdcKalTrackCol();
      DstMdcKalTrackCol::iterator iter_mdc       = dstMdcKalTrackCol->begin();
      trackID                                    = 0;
      for ( ; iter_mdc != dstMdcKalTrackCol->end(); iter_mdc++, trackID++ )
      {
        RecMdcKalTrack* recMdcKalTrack = new RecMdcKalTrack();
        *recMdcKalTrack                = **iter_mdc;
        ( *mdcKalTrackCol ).push_back( recMdcKalTrack );
        log << MSG::INFO << " MdcKalTrack ID = " << trackID
            << " MdcKalTrack Nster = " << ( *iter_mdc )->nster()
            << " MdcKalTrack poca = " << ( *iter_mdc )->poca() << endmsg;
      }
      sc = dataSvc->registerObject( EventModel::Recon::RecMdcKalTrackCol, mdcKalTrackCol );
    }
  }
 
  //**********************************************************************
  // now convert
  TIter        dedxIter( m_recMdcDedxCol );
  TRecMdcDedx* recMdcDedxRoot = 0;
  while ( ( recMdcDedxRoot = (TRecMdcDedx*)dedxIter.Next() ) )
  {
    double dedxHit    = recMdcDedxRoot->dedxHit();
    double dedxEsat   = recMdcDedxRoot->dedxEsat();
    double dedxNoRun  = recMdcDedxRoot->dedxNoRun();
    double dedxMoment = recMdcDedxRoot->dedxMoment();
 
    int trackId    = recMdcDedxRoot->trackId();
    int particleId = recMdcDedxRoot->particleId();
    int status     = recMdcDedxRoot->status();
    int truncAlg   = recMdcDedxRoot->truncAlg();
    // double pb[5];
    // for ( int i = 0; i < 5; i++)
    //      pb[i] = dedxRoot->prob(i);
 
    // double numSigmaE  = dedxRoot->numSigmaE();
    // double numSigmaMu = dedxRoot->numSigmaMu();
    // double numSigmaPi = dedxRoot->numSigmaPi();
    // double numSigmaK = dedxRoot->numSigmaK();
    // double numSigmaP = dedxRoot->numSigmaP();
    double chi[5];
    chi[0]           = recMdcDedxRoot->chiE();
    chi[1]           = recMdcDedxRoot->chiMu();
    chi[2]           = recMdcDedxRoot->chiPi();
    chi[3]           = recMdcDedxRoot->chiK();
    chi[4]           = recMdcDedxRoot->chiP();
    int numGoodHits  = recMdcDedxRoot->numGoodHits();
    int numTotalHits = recMdcDedxRoot->numTotalHits();
 
    double probPH   = recMdcDedxRoot->probPH();
    double normPH   = recMdcDedxRoot->normPH();
    double errorPH  = recMdcDedxRoot->errorPH();
    double twentyPH = recMdcDedxRoot->twentyPH();
    double dedxExpect[5], sigmaDedx[5], pidProb[5];
 
    for ( int i = 0; i < 5; i++ )
    {
      dedxExpect[i] = recMdcDedxRoot->dedxExpect( i );
      sigmaDedx[i]  = recMdcDedxRoot->sigmaDedx( i );
      pidProb[i]    = recMdcDedxRoot->pidProb( i );
    }
 
    log << MSG::DEBUG << "TObjectToDataObject: check Reconstrunction of dE/dx root::"
        << " trackId: " << trackId << " particleId: " << particleId << " status: " << status
        << " truncAlg: " << truncAlg << " chi[2]: " << chi[2]
        << " numTotalHits: " << numTotalHits << " probPH: " << probPH
        << " errorPH: " << errorPH << " dedxExpect[2]: " << dedxExpect[2] << endmsg;
 
    RecMdcDedx* recMdcDedx                       = new RecMdcDedx();
    m_common.m_rootRecMdcDedxMap[recMdcDedxRoot] = recMdcDedx;
 
    recMdcDedx->setTrackId( trackId );
    recMdcDedx->setParticleId( particleId );
    recMdcDedx->setStatus( status );
    recMdcDedx->setTruncAlg( truncAlg );
    // dedxTds->setProb(pb);
    // dedxTds->setNumSigmaE(numSigmaE);
    // dedxTds->setNumSigmaMu(numSigmaMu);
    // dedxTds->setNumSigmaPi(numSigmaPi);
    // dedxTds->setNumSigmaK(numSigmaK);
    // dedxTds->setNumSigmaP(numSigmaP);
    recMdcDedx->setChi( chi );
    recMdcDedx->setNumGoodHits( numGoodHits );
    recMdcDedx->setNumTotalHits( numTotalHits );
 
    recMdcDedx->setProbPH( probPH );
    recMdcDedx->setNormPH( normPH );
    recMdcDedx->setErrorPH( errorPH );
    recMdcDedx->setTwentyPH( twentyPH );
    // for (int i=0; i<5; i++){
    recMdcDedx->setDedxExpect( dedxExpect );
    recMdcDedx->setSigmaDedx( sigmaDedx );
    recMdcDedx->setPidProb( pidProb );
 
    recMdcDedx->setDedxHit( dedxHit );
    recMdcDedx->setDedxEsat( dedxEsat );
    recMdcDedx->setDedxNoRun( dedxNoRun );
    recMdcDedx->setDedxMoment( dedxMoment );
    //}
 
    DedxHitRefVec              theDedxHitRefVec;
    RecMdcDedxHitCol::iterator iter = recMdcDedxHitCol->begin();
    for ( ; iter != recMdcDedxHitCol->end(); iter++ )
    {
      if ( ( *iter )->getTrkId() == trackId )
      {
        SmartRef<RecMdcDedxHit> refDedxHit( *iter );
        theDedxHitRefVec.push_back( refDedxHit );
      }
    }
    recMdcDedx->setVecDedxHits( theDedxHitRefVec );
    int nhits = recMdcDedx->getVecDedxHits().size();
    // std::cout<<" mdc hits: "<<nhits<< std::endl;
    for ( int ii = 0; ii < nhits; ii++ )
    {
      Identifier id( recMdcDedx->getVecDedxHits()[ii]->getMdcId() );
      int        layer = MdcID::layer( id );
      int        wire  = MdcID::wire( id );
      // cout<<"layer: "<<layer<<" wire: "<<wire<<endl;
    }
 
    int mdcTrackId = recMdcDedxRoot->mdcTrackId();
    // std::cout << __FILE__ << __LINE__ << " size: " << mdcTrackCol->size() << "  id: " <<
    // mdcTrackId << std::endl;
    if ( mdcTrackId >= 0 )
    {
      recMdcDedx->setMdcTrack(
          dynamic_cast<RecMdcTrack*>( mdcTrackCol->containedObject( mdcTrackId ) ) );
    }
    int mdcKalTrackId = recMdcDedxRoot->mdcKalTrackId();
    // std::cout << __FILE__ << __LINE__ << " size: " << mdcKalTrackCol->size() << "  id: " <<
    // mdcKalTrackId<< std::endl;
    if ( mdcKalTrackId >= 0 )
    {
      recMdcDedx->setMdcKalTrack(
          dynamic_cast<RecMdcKalTrack*>( mdcKalTrackCol->containedObject( mdcKalTrackId ) ) );
    }
 
    recMdcDedxCol->push_back( recMdcDedx );
    // delete dedxTds;
    //  dedxTds = NULL;
  }
  // m_dedxCol->Delete();
  delete m_recMdcDedxCol;
  m_recMdcDedxCol = 0;
 
  return StatusCode::SUCCESS;
}
 
StatusCode RecMdcDedxCnv::DataObjectToTObject( DataObject* obj, RootAddress* rootaddr ) {
  MsgStream log( msgSvc(), "RecMdcDedxCnv" );
  log << MSG::DEBUG << "RecMdcDedxCnv::DataObjectToTObject" << endmsg;
  StatusCode sc = StatusCode::SUCCESS;
 
  RecMdcDedxCol* recMdcDedxCol = dynamic_cast<RecMdcDedxCol*>( obj );
  if ( !recMdcDedxCol )
  {
    log << MSG::ERROR << "Could not downcast to RecMdcDedxCol" << endmsg;
    return StatusCode::FAILURE;
  }
 
  DataObject* evt;
  m_eds->findObject( EventModel::Recon::Event, evt );
  if ( evt == NULL )
  {
    log << MSG::ERROR << "Could not get RecEvent in TDS " << endmsg;
    return StatusCode::FAILURE;
  }
  ReconEvent* devtTds = dynamic_cast<ReconEvent*>( evt );
  if ( !devtTds )
  { log << MSG::ERROR << "RecMdcDedxCnv:Could not downcast to TDS RecEvent" << endmsg; }
 
  IOpaqueAddress* addr;
  m_cnvSvc->getRecTrackCnv()->createRep( evt, addr );
  TRecTrackEvent* recEvt = m_cnvSvc->getRecTrackCnv()->getWriteObject();
 
  const TObjArray* m_recMdcDedxCol = recEvt->getRecMdcDedxCol();
  if ( !m_recMdcDedxCol ) return sc;
  recEvt->clearRecMdcDedxCol(); // necessary in case there is I/O at the same time since array
                                // is static
 
  //******************** get according mdcTrack and MdcKalmanTrack of MdcDedx  ****************
 
  RecMdcTrackCol::iterator    recMdcTrackColbegin, recMdcTrackColend;
  RecMdcKalTrackCol::iterator recMdcKalTrackColbegin, recMdcKalTrackColend;
 
  IDataProviderSvc* dataSvc = 0;
  sc = serviceLocator()->getService( "EventDataSvc", IDataProviderSvc::interfaceID(),
                                     (IInterface*&)dataSvc );
  if ( sc.isFailure() )
  {
    log << MSG::FATAL << "Could not get EventDataSvc in RecMdcDedxCnv" << endmsg;
    return sc;
  }
  SmartDataPtr<RecMdcTrackCol> recMdcTrackCol( dataSvc, EventModel::Recon::RecMdcTrackCol );
  if ( !recMdcTrackCol )
  {
    log << MSG::ERROR << "Could not downcast to RecMdcTrackCol" << endmsg;
    return StatusCode::FAILURE;
  }
  else
  {
    recMdcTrackColbegin = recMdcTrackCol->begin();
    recMdcTrackColend   = recMdcTrackCol->end();
  }
 
  SmartDataPtr<RecMdcKalTrackCol> recMdcKalTrackCol( dataSvc,
                                                     EventModel::Recon::RecMdcKalTrackCol );
  if ( !recMdcKalTrackCol )
  {
    log << MSG::ERROR << "Could not downcast to RecMdcKalTrackCol" << endmsg;
    return StatusCode::FAILURE;
  }
  else
  {
    recMdcKalTrackColbegin = recMdcKalTrackCol->begin();
    recMdcKalTrackColend   = recMdcKalTrackCol->end();
  }
  //******************************************************************************
 
  RecMdcDedxCol::const_iterator recMdcDedx;
 
  for ( recMdcDedx = recMdcDedxCol->begin(); recMdcDedx != recMdcDedxCol->end(); recMdcDedx++ )
  {
    Int_t trackId    = ( *recMdcDedx )->trackId();
    Int_t particleId = ( *recMdcDedx )->particleId();
    Int_t status     = ( *recMdcDedx )->status();
    Int_t truncAlg   = ( *recMdcDedx )->truncAlg();
    // Double_t pb[5];
    // for (Int_t i = 0; i < 5; i++)
    // pb[i] = (*dedxTds)->prob(i);
    Double_t dedxHit    = ( *recMdcDedx )->getDedxHit();
    Double_t dedxEsat   = ( *recMdcDedx )->getDedxEsat();
    Double_t dedxNoRun  = ( *recMdcDedx )->getDedxNoRun();
    Double_t dedxMoment = ( *recMdcDedx )->getDedxMoment();
 
    Double_t chiE  = ( *recMdcDedx )->chi( 0 );
    Double_t chiMu = ( *recMdcDedx )->chi( 1 );
    Double_t chiPi = ( *recMdcDedx )->chi( 2 );
    Double_t chiK  = ( *recMdcDedx )->chi( 3 );
    Double_t chiP  = ( *recMdcDedx )->chi( 4 );
 
    Int_t numGoodHits  = ( *recMdcDedx )->numGoodHits();
    Int_t numTotalHits = ( *recMdcDedx )->numTotalHits();
 
    Double_t probPH   = ( *recMdcDedx )->probPH();
    Double_t normPH   = ( *recMdcDedx )->normPH();
    Double_t errorPH  = ( *recMdcDedx )->errorPH();
    Double_t twentyPH = ( *recMdcDedx )->twentyPH();
    // Double_t fracErrPH = (*dedxTds)-> fracErrPH();
    // Double_t minIronPH = (*dedxTds)->minIronPH();
    // Double_t corrPH = (*dedxTds)->corrPH();
    double dedxExpect[5], sigmaDedx[5], pidProb[5], chi[5];
    for ( int i = 0; i < 5; i++ )
    {
      chi[i]        = ( *recMdcDedx )->chi( i );
      dedxExpect[i] = ( *recMdcDedx )->getDedxExpect( i );
      sigmaDedx[i]  = ( *recMdcDedx )->getSigmaDedx( i );
      pidProb[i]    = ( *recMdcDedx )->getPidProb( i );
    }
 
    log << MSG::DEBUG << "DataObjectToTObject: check Reconstrunction of dE/dx::"
        << " trackId: " << trackId << " particleId: " << particleId << " status: " << status
        << " truncAlg: " << truncAlg << " chiPi: " << chiPi
        << " numTotalHits: " << numTotalHits << " probPH: " << probPH
        << " errorPH: " << errorPH << " dedxExpect[2]: " << dedxExpect[2] << endmsg;
 
    TRecMdcDedx* recMdcDedxRoot = new TRecMdcDedx();
    // m_common.m_recMdcDedxMap[(*recMdcDedx)] = recMdcDedxRoot;
    // std::cout<<"check write to Reconstrunction root particle Id is
    // "<<dedxRoot->particleId()<<endl;
    recMdcDedxRoot->setTrackId( trackId );
    recMdcDedxRoot->setParticleId( particleId );
    recMdcDedxRoot->setStatus( status );
    recMdcDedxRoot->setTruncAlg( truncAlg );
 
    // dedxRoot->setProb(pb);
    recMdcDedxRoot->setChiE( chiE );
    recMdcDedxRoot->setChiMu( chiMu );
    recMdcDedxRoot->setChiPi( chiPi );
    recMdcDedxRoot->setChiK( chiK );
    recMdcDedxRoot->setChiP( chiP );
 
    recMdcDedxRoot->setNumGoodHits( numGoodHits );
    recMdcDedxRoot->setNumTotalHits( numTotalHits );
 
    recMdcDedxRoot->setProbPH( probPH );
    recMdcDedxRoot->setNormPH( normPH );
    recMdcDedxRoot->setErrorPH( errorPH );
    recMdcDedxRoot->setTwentyPH( twentyPH );
    //  for (int i=0; i<5; i++){
    recMdcDedxRoot->setChi( chi );
    recMdcDedxRoot->setDedxExpect( dedxExpect );
    recMdcDedxRoot->setSigmaDedx( sigmaDedx );
    recMdcDedxRoot->setPidProb( pidProb );
 
    recMdcDedxRoot->setDedxHit( dedxHit );
    recMdcDedxRoot->setDedxEsat( dedxEsat );
    recMdcDedxRoot->setDedxNoRun( dedxNoRun );
    recMdcDedxRoot->setDedxMoment( dedxMoment );
    //}
 
    if ( ( *recMdcDedx )->isMdcTrackValid() )
    {
      RecMdcTrackCol::iterator it =
          find( recMdcTrackColbegin, recMdcTrackColend, ( *recMdcDedx )->getMdcTrack() );
      recMdcDedxRoot->setMdcTrackId( it - recMdcTrackColbegin );
    }
    else { recMdcDedxRoot->setMdcTrackId( -1 ); }
 
    if ( ( *recMdcDedx )->isMdcKalTrackValid() )
    {
      RecMdcKalTrackCol::iterator it = find( recMdcKalTrackColbegin, recMdcKalTrackColend,
                                             ( *recMdcDedx )->getMdcKalTrack() );
      recMdcDedxRoot->setMdcKalTrackId( it - recMdcKalTrackColbegin );
    }
    else { recMdcDedxRoot->setMdcKalTrackId( -1 ); }
 
    log << MSG::INFO << "check Reconstrunction root"
        << " particle Id is : " << particleId << "track id is : " << trackId
        << " and status is : " << status << endmsg;
 
    // dedxRoot->setFracErrPH(fracErrPH);
    // dedxRoot->setMinIronPH(minIronPH);
    // dedxRoot->setCorrPH(corrPH);
    recEvt->addRecMdcDedx( recMdcDedxRoot );
  }
 
  return StatusCode::SUCCESS;
}
 
#endif