Single.cxx

Go to the documentation of this file.

#include "GaudiKernel/Bootstrap.h"
#include "GaudiKernel/IDataProviderSvc.h"
#include "GaudiKernel/ISvcLocator.h"
#include "GaudiKernel/MsgStream.h"
#include "GaudiKernel/PropertyMgr.h"
#include "GaudiKernel/SmartDataPtr.h"
#include "VertexDbSvc/IVertexDbSvc.h"
 
#include "EventModel/Event.h"
#include "EventModel/EventModel.h"
 
#include "DstEvent/TofHitStatus.h"
#include "EventModel/EventHeader.h"
#include "EvtRecEvent/EvtRecEvent.h"
#include "EvtRecEvent/EvtRecTrack.h"
#include "McTruth/McParticle.h"
#include "MucRecEvent/RecMucTrack.h"
 
#include "CLHEP/Vector/LorentzVector.h"
#include "CLHEP/Vector/ThreeVector.h"
#include "CLHEP/Vector/TwoVector.h"
#include "GaudiKernel/Bootstrap.h"
#include "GaudiKernel/IHistogramSvc.h"
#include "GaudiKernel/INTupleSvc.h"
#include "GaudiKernel/NTuple.h"
#include "TMath.h"
using CLHEP::Hep2Vector;
using CLHEP::Hep3Vector;
using CLHEP::HepLorentzVector;
#include "CLHEP/Geometry/Point3D.h"
#ifndef ENABLE_BACKWARDS_COMPATIBILITY
typedef HepGeom::Point3D<double> HepPoint3D;
#endif
#include "Single.h"
 
#include "ParticleID/ParticleID.h"
#include "VertexFit/Helix.h"
#include "VertexFit/KinematicFit.h"
#include "VertexFit/VertexFit.h"
 
#include <vector>
 
typedef std::vector<int>              Vint;
typedef std::vector<HepLorentzVector> Vp4;
 
int Ncut0, Ncut1, Ncut2, Ncut3, Ncut4;
 
/////////////////////////////////////////////////////////////////////////////
DECLARE_COMPONENT( Single )
Single::Single( const std::string& name, ISvcLocator* pSvcLocator )
    : Algorithm( name, pSvcLocator ) {
  // Declare the properties
  declareProperty( "Vxy0cut", m_vxy0cut = 2.0 );
  declareProperty( "Vz0cut", m_vz0cut = 20.0 );
  declareProperty( "angdcut", m_angdcut = 30.0 );
  declareProperty( "CheckDedx", m_checkDedx = 1 );
  declareProperty( "CheckTof", m_checkTof = 1 );
  declareProperty( "CheckExt", m_checkExt = 1 );
  declareProperty( "CheckEmc", m_checkEmc = 1 );
  declareProperty( "CheckMuc", m_checkMuc = 1 );
  declareProperty( "ParticleID", m_particleID = -211 );
  declareProperty( "Itrack", m_itrack = 0 );
  declareProperty( "PID_mode", m_pid_mode = 0 ); // 0:Combine  1:Only Dedx  2:Only Tof
}
 
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
StatusCode Single::initialize() {
  MsgStream log( msgSvc(), name() );
 
  log << MSG::INFO << "in initialize()" << endmsg;
 
  StatusCode status;
  NTuplePtr  nt1( ntupleSvc(), "FILE1/mdctrack" );
  if ( nt1 ) m_tuple1 = nt1;
  else
  {
    m_tuple1 =
        ntupleSvc()->book( "FILE1/mdctrack", CLID_ColumnWiseTuple, "ks N-Tuple example" );
    if ( m_tuple1 )
    {
      status = m_tuple1->addItem( "mdc_vx0", m_vx0 );
      status = m_tuple1->addItem( "mdc_vy0", m_vy0 );
      status = m_tuple1->addItem( "mdc_vz0", m_vz0 );
      status = m_tuple1->addItem( "mdc_pch", m_pch );
      status = m_tuple1->addItem( "mdc_trackid", m_trackid );
      status = m_tuple1->addItem( "mdc_charge", m_charge );
      status = m_tuple1->addItem( "mdc_pxy0", m_pxy0 );
      status = m_tuple1->addItem( "mdc_px0", m_px0 );
      status = m_tuple1->addItem( "mdc_py0", m_py0 );
      status = m_tuple1->addItem( "mdc_pz0", m_pz0 );
      status = m_tuple1->addItem( "mdc_r0", m_r0 );
      status = m_tuple1->addItem( "mdc_phi", m_phi );
      status = m_tuple1->addItem( "mdc_theta", m_theta );
      status = m_tuple1->addItem( "mdc_ndof", m_ndof );
      status = m_tuple1->addItem( "mdc_nster", m_nster );
      status = m_tuple1->addItem( "mdc_stat", m_stat );
      status = m_tuple1->addItem( "mdc_angd", m_angd );
      status = m_tuple1->addItem( "mdc_rvxy0", m_rvxy0 );
      status = m_tuple1->addItem( "mdc_rvz0", m_rvz0 );
      status = m_tuple1->addItem( "mdc_d0", m_d0 );
      status = m_tuple1->addItem( "mdc_phi0", m_phi0 );
      status = m_tuple1->addItem( "mdc_kappa", m_kappa );
      status = m_tuple1->addItem( "mdc_dz", m_dzhelix );
      status = m_tuple1->addItem( "mdc_tanlamda", m_tanlamda );
      status = m_tuple1->addItem( "mdc_err11", m_err11 );
      status = m_tuple1->addItem( "mdc_err21", m_err21 );
      status = m_tuple1->addItem( "mdc_err22", m_err22 );
      status = m_tuple1->addItem( "mdc_err31", m_err31 );
      status = m_tuple1->addItem( "mdc_err32", m_err32 );
      status = m_tuple1->addItem( "mdc_err33", m_err33 );
      status = m_tuple1->addItem( "mdc_err41", m_err41 );
      status = m_tuple1->addItem( "mdc_err42", m_err42 );
      status = m_tuple1->addItem( "mdc_err43", m_err43 );
      status = m_tuple1->addItem( "mdc_err44", m_err44 );
      status = m_tuple1->addItem( "mdc_err51", m_err51 );
      status = m_tuple1->addItem( "mdc_err52", m_err52 );
      status = m_tuple1->addItem( "mdc_err53", m_err53 );
      status = m_tuple1->addItem( "mdc_err54", m_err54 );
      status = m_tuple1->addItem( "mdc_err55", m_err55 );
 
      status = m_tuple1->addItem( "mdc_kal_px", m_mdc_kal_px );
      status = m_tuple1->addItem( "mdc_kal_py", m_mdc_kal_py );
      status = m_tuple1->addItem( "mdc_kal_pz", m_mdc_kal_pz );
      status = m_tuple1->addItem( "mdc_kal_p", m_mdc_kal_p );
      status = m_tuple1->addItem( "mdc_kal_pxy", m_mdc_kal_pxy );
      status = m_tuple1->addItem( "mdc_kal_costheta", m_mdc_kal_costheta );
 
      status = m_tuple1->addItem( "NGch", m_ngch );
    }
    else
    {
      log << MSG::ERROR << "    Cannot book N-tuple:" << long( m_tuple1 ) << endmsg;
      return StatusCode::FAILURE;
    }
  }
 
  // check mdctrack
 
  if ( m_checkDedx == 1 )
  {
    NTuplePtr nt2( ntupleSvc(), "FILE1/dedx" );
    if ( nt2 ) m_tuple2 = nt2;
    else
    {
      m_tuple2 = ntupleSvc()->book( "FILE1/dedx", CLID_ColumnWiseTuple, "ks N-Tuple example" );
      if ( m_tuple2 )
      {
        status = m_tuple2->addItem( "dedx_ptrk", m_ptrk );
        status = m_tuple2->addItem( "dedx_chie", m_chie );
        status = m_tuple2->addItem( "dedx_chimu", m_chimu );
        status = m_tuple2->addItem( "dedx_chipi", m_chipi );
        status = m_tuple2->addItem( "dedx_chik", m_chik );
        status = m_tuple2->addItem( "dedx_chip", m_chip );
        status = m_tuple2->addItem( "dedx_probPH", m_probPH );
        status = m_tuple2->addItem( "dedx_normPH", m_normPH );
        status = m_tuple2->addItem( "dedx_ghit", m_ghit );
        status = m_tuple2->addItem( "dedx_thit", m_thit );
      }
      else
      {
        log << MSG::ERROR << "    Cannot book N-tuple:" << long( m_tuple2 ) << endmsg;
        return StatusCode::FAILURE;
      }
    }
  }
  // check dE/dx
 
  if ( m_checkExt == 1 )
  {
    NTuplePtr nt3( ntupleSvc(), "FILE1/ext" );
    if ( nt3 ) m_tuple3 = nt3;
    else
    {
      m_tuple3 = ntupleSvc()->book( "FILE1/ext", CLID_ColumnWiseTuple, "ks N-Tuple example" );
      if ( m_tuple3 )
      {
 
        status = m_tuple3->addItem( "ext_tof1", m_tof1 );
        status = m_tuple3->addItem( "ext_tof1path", m_tof1path );
        status = m_tuple3->addItem( "ext_tof1z", m_tof1z );
        status = m_tuple3->addItem( "ext_tof1t", m_tof1t );
        status = m_tuple3->addItem( "ext_tof1x", m_tof1x );
        status = m_tuple3->addItem( "ext_tof1y", m_tof1y );
 
        status = m_tuple3->addItem( "ext_tof2", m_tof2 );
        status = m_tuple3->addItem( "ext_tof2path", m_tof2path );
        status = m_tuple3->addItem( "ext_tof2z", m_tof2z );
        status = m_tuple3->addItem( "ext_tof2t", m_tof2t );
        status = m_tuple3->addItem( "ext_tof2x", m_tof2x );
        status = m_tuple3->addItem( "ext_tof2y", m_tof2y );
 
        status = m_tuple3->addItem( "ext_emctheta", m_emctheta );
        status = m_tuple3->addItem( "ext_emcphi", m_emcphi );
        status = m_tuple3->addItem( "ext_emcpath", m_emcpath );
 
        status = m_tuple3->addItem( "ext_mucz", m_mucz );
        status = m_tuple3->addItem( "ext_muct", m_muct );
        status = m_tuple3->addItem( "ext_mucx", m_mucx );
        status = m_tuple3->addItem( "ext_mucy", m_mucy );
      }
      else
      {
        log << MSG::ERROR << "    Cannot book N-tuple:" << long( m_tuple3 ) << endmsg;
        return StatusCode::FAILURE;
      }
    }
  }
  // check exttrack
 
  if ( m_checkTof == 1 )
  {
    NTuplePtr nt4( ntupleSvc(), "FILE1/tof" );
    if ( nt4 ) m_tuple4 = nt4;
    else
    {
      m_tuple4 = ntupleSvc()->book( "FILE1/tof", CLID_ColumnWiseTuple, "ks N-Tuple example" );
      if ( m_tuple2 )
      {
        status = m_tuple4->addItem( "tof_path", m_path );
        status = m_tuple4->addItem( "tof_zrhit", m_zrhit );
        status = m_tuple4->addItem( "tof_ph", m_ph );
        status = m_tuple4->addItem( "tof_tof", m_tof );
        status = m_tuple4->addItem( "tof_errtof", m_errtof );
        status = m_tuple4->addItem( "tof_beta", m_beta );
        status = m_tuple4->addItem( "tof_texpe", m_texpe );
        status = m_tuple4->addItem( "tof_texpmu", m_texpmu );
        status = m_tuple4->addItem( "tof_texppi", m_texppi );
        status = m_tuple4->addItem( "tof_texpka", m_texpka );
        status = m_tuple4->addItem( "tof_texppro", m_texppro );
        status = m_tuple4->addItem( "tof_toffsete", m_toffsete );
        status = m_tuple4->addItem( "tof_toffsetmu", m_toffsetmu );
        status = m_tuple4->addItem( "tof_toffsetpi", m_toffsetpi );
        status = m_tuple4->addItem( "tof_toffsetka", m_toffsetka );
        status = m_tuple4->addItem( "tof_toffsetpro", m_toffsetpro );
        status = m_tuple4->addItem( "tof_toffsetatpr", m_toffsetatpr );
        status = m_tuple4->addItem( "tof_sigmae", m_sigmae );
        status = m_tuple4->addItem( "tof_sigmamu", m_sigmamu );
        status = m_tuple4->addItem( "tof_sigmapi", m_sigmapi );
        status = m_tuple4->addItem( "tof_sigmaka", m_sigmaka );
        status = m_tuple4->addItem( "tof_sigmapro", m_sigmapro );
        status = m_tuple4->addItem( "tof_sigmaatpr", m_sigmaatpr );
        status = m_tuple4->addItem( "tof_quality", m_quality );
        status = m_tuple4->addItem( "tof_t0", m_t0 );
        status = m_tuple4->addItem( "tof_errt0", m_errt0 );
        status = m_tuple4->addItem( "tof_errz", m_errz );
        status = m_tuple4->addItem( "tof_phi", m_phitof );
        status = m_tuple4->addItem( "tof_errphi", m_errphi );
        status = m_tuple4->addItem( "tof_energy", m_energy );
        status = m_tuple4->addItem( "tof_errenergy", m_errenergy );
      }
      else
      {
        log << MSG::ERROR << "    Cannot book N-tuple:" << long( m_tuple4 ) << endmsg;
        return StatusCode::FAILURE;
      }
    }
  }
  // check tof
 
  if ( m_checkEmc == 1 )
  {
    NTuplePtr nt5( ntupleSvc(), "FILE1/emc" );
    if ( nt5 ) m_tuple5 = nt5;
    else
    {
      m_tuple5 = ntupleSvc()->book( "FILE1/emc", CLID_ColumnWiseTuple, "ks N-Tuple example" );
      if ( m_tuple5 )
      {
 
        status = m_tuple5->addItem( "emc_x", m_x );
        status = m_tuple5->addItem( "emc_y", m_y );
        status = m_tuple5->addItem( "emc_z", m_z );
        status = m_tuple5->addItem( "emc_theta", m_thetaemc );
        status = m_tuple5->addItem( "emc_phi", m_phiemc );
        status = m_tuple5->addItem( "emc_dx", m_dx );
        status = m_tuple5->addItem( "emc_dy", m_dy );
        status = m_tuple5->addItem( "emc_dz", m_dz );
        status = m_tuple5->addItem( "emc_dtheta", m_dtheta );
        status = m_tuple5->addItem( "emc_dphi", m_dphi );
        status = m_tuple5->addItem( "emc_energy", m_energyemc );
        status = m_tuple5->addItem( "emc_dE", m_de );
        status = m_tuple5->addItem( "emc_eseed", m_eseed );
        status = m_tuple5->addItem( "emc_e3x3", m_e3x3 );
        status = m_tuple5->addItem( "emc_e5x5", m_e5x5 );
        status = m_tuple5->addItem( "emc_secp", m_secp );
        status = m_tuple5->addItem( "emc_latp", m_latp );
      }
      else
      {
        log << MSG::ERROR << "    Cannot book N-tuple:" << long( m_tuple5 ) << endmsg;
        return StatusCode::FAILURE;
      }
    }
  }
  // check emc
 
  if ( m_checkMuc == 1 )
  {
    NTuplePtr nt6( ntupleSvc(), "FILE1/muc" );
    if ( nt6 ) m_tuple6 = nt6;
    else
    {
      m_tuple6 = ntupleSvc()->book( "FILE1/muc", CLID_ColumnWiseTuple, "ks N-Tuple example" );
      if ( m_tuple6 )
      {
 
        status = m_tuple6->addItem( "muc_depth", m_depth );
        status = m_tuple6->addItem( "muc_chi2", m_chi2 );
        status = m_tuple6->addItem( "muc_rms", m_rms );
        status = m_tuple6->addItem( "muc_xpos", m_xpos );
        status = m_tuple6->addItem( "muc_ypos", m_ypos );
        status = m_tuple6->addItem( "muc_zpos", m_zpos );
        status = m_tuple6->addItem( "muc_xsigma", m_xpossigma );
        status = m_tuple6->addItem( "muc_ysigma", m_ypossigma );
        status = m_tuple6->addItem( "muc_zsigma", m_zpossigma );
        status = m_tuple6->addItem( "muc_px", m_px );
        status = m_tuple6->addItem( "muc_py", m_py );
        status = m_tuple6->addItem( "muc_pz", m_pz );
        status = m_tuple6->addItem( "muc_distance", m_distance );
        status = m_tuple6->addItem( "muc_deltaphi", m_deltaphi );
      }
      else
      {
        log << MSG::ERROR << "    Cannot book N-tuple:" << long( m_tuple6 ) << endmsg;
        return StatusCode::FAILURE;
      }
    }
  }
  // check muc
 
  NTuplePtr nt7( ntupleSvc(), "FILE1/runinfo" );
  if ( nt7 ) m_tuple7 = nt7;
  else
  {
    m_tuple7 =
        ntupleSvc()->book( "FILE1/runinfo", CLID_ColumnWiseTuple, "ks N-Tuple example" );
    if ( m_tuple7 )
    {
 
      status = m_tuple7->addItem( "runinfo_runNo", m_runNo );
      status = m_tuple7->addItem( "runinfo_event", m_event );
      status = m_tuple7->addItem( "runinfo_totcharged", m_totcharged );
      status = m_tuple7->addItem( "runinfo_totneutral", m_totneutral );
      status = m_tuple7->addItem( "runinfo_tottracks", m_tottracks );
      status = m_tuple7->addItem( "runinfo_ngood", m_ngood );
      status = m_tuple7->addItem( "runinfo_ncharge", m_ncharge );
      status = m_tuple7->addItem( "kal_n", m_kaln );
      status = m_tuple7->addItem( "costheta_n", m_costhetan );
      status = m_tuple7->addItem( "p_n", m_pn );
      status = m_tuple7->addItem( "charge", m_mat_charge );
      status = m_tuple7->addItem( "diff_pxy", m_diff_pxy );
      status = m_tuple7->addItem( "diff_p", m_diff_pch );
      status = m_tuple7->addItem( "diff_costheta", m_diff_costheta );
    }
    else
    {
      log << MSG::ERROR << "    Cannot book N-tuple:" << long( m_tuple7 ) << endmsg;
      return StatusCode::FAILURE;
    }
  }
  // check runinfo
 
  NTuplePtr nt8( ntupleSvc(), "FILE1/mcpart" );
  if ( nt8 ) m_tuple8 = nt8;
  else
  {
    m_tuple8 = ntupleSvc()->book( "FILE1/mcpart", CLID_ColumnWiseTuple, "ks N-Tuple example" );
    if ( m_tuple8 )
    {
      status = m_tuple8->addItem( "mcpart_vx", m_vx );
      status = m_tuple8->addItem( "mcpart_vy", m_vy );
      status = m_tuple8->addItem( "mcpart_vz", m_vz );
      status = m_tuple8->addItem( "mcpart_vr0", m_vr0 );
      status = m_tuple8->addItem( "mcpart_px", m_pxmc );
      status = m_tuple8->addItem( "mcpart_py", m_pymc );
      status = m_tuple8->addItem( "mcpart_pz", m_pzmc );
      status = m_tuple8->addItem( "mcpart_e", m_e );
      status = m_tuple8->addItem( "mcpart_p", m_p );
      status = m_tuple8->addItem( "mcpart_pxy", m_pxymc );
      status = m_tuple8->addItem( "mcpart_theta", m_thetamc );
      status = m_tuple8->addItem( "mcpart_costheta", m_costhetamc );
 
      status = m_tuple8->addItem( "mcpart_phi", m_phimc );
      status = m_tuple8->addItem( "mcpart_m", m_m );
      status = m_tuple8->addItem( "mcpart_pro", m_pro );
      status = m_tuple8->addItem( "mcpart_index", m_index );
      status = m_tuple8->addItem( "totalcharged", m_mctotcharged );
    }
    else
    {
      log << MSG::ERROR << "    Cannot book N-tuple:" << long( m_tuple8 ) << endmsg;
      return StatusCode::FAILURE;
    }
  }
  // check mctruth
 
  NTuplePtr nt9( ntupleSvc(), "FILE1/pid_kal" );
  if ( nt9 ) m_tuple9 = nt9;
  else
  {
    m_tuple9 =
        ntupleSvc()->book( "FILE1/pid_kal", CLID_ColumnWiseTuple, "ks N-Tuple example" );
    if ( m_tuple9 )
    {
      status = m_tuple9->addItem( "pid_ptrk", m_ptrk_pid );
      status = m_tuple9->addItem( "pid_cost", m_cost_pid );
      status = m_tuple9->addItem( "pid_dedx", m_dedx_pid );
      status = m_tuple9->addItem( "pid_tof1", m_tof1_pid );
      status = m_tuple9->addItem( "pid_tof2", m_tof2_pid );
      status = m_tuple9->addItem( "prob_pion", m_prob_pid_pion );
      status = m_tuple9->addItem( "prob_kaon", m_prob_pid_kaon );
      status = m_tuple9->addItem( "prob_proton", m_prob_pid_proton );
 
      status = m_tuple9->addItem( "kal_px", m_kalpx );
      status = m_tuple9->addItem( "kal_py", m_kalpy );
      status = m_tuple9->addItem( "kal_pz", m_kalpz );
      status = m_tuple9->addItem( "kal_p", m_kalp );
      status = m_tuple9->addItem( "kal_pxy", m_kalpxy );
 
      status = m_tuple9->addItem( "kal_costheta", m_kalcostheta );
 
      status = m_tuple9->addItem( "kal_d0", m_d0kal );
      status = m_tuple9->addItem( "kal_phi0", m_phi0kal );
      status = m_tuple9->addItem( "kal_kappa", m_kappakal );
      status = m_tuple9->addItem( "kal_dz", m_dzhelixkal );
      status = m_tuple9->addItem( "kal_tanlamda", m_tanlamdakal );
      status = m_tuple9->addItem( "kal_err11", m_err11kal );
      status = m_tuple9->addItem( "kal_err21", m_err21kal );
      status = m_tuple9->addItem( "kal_err22", m_err22kal );
      status = m_tuple9->addItem( "kal_err31", m_err31kal );
      status = m_tuple9->addItem( "kal_err32", m_err32kal );
      status = m_tuple9->addItem( "kal_err33", m_err33kal );
      status = m_tuple9->addItem( "kal_err41", m_err41kal );
      status = m_tuple9->addItem( "kal_err42", m_err42kal );
      status = m_tuple9->addItem( "kal_err43", m_err43kal );
      status = m_tuple9->addItem( "kal_err44", m_err44kal );
      status = m_tuple9->addItem( "kal_err51", m_err51kal );
      status = m_tuple9->addItem( "kal_err52", m_err52kal );
      status = m_tuple9->addItem( "kal_err53", m_err53kal );
      status = m_tuple9->addItem( "kal_err54", m_err54kal );
      status = m_tuple9->addItem( "kal_err55", m_err55kal );
      status = m_tuple9->addItem( "kal_nn", m_kalnn );
      status = m_tuple9->addItem( "costheta_nn", m_costhetann );
      status = m_tuple9->addItem( "p_nn", m_pnn );
 
      status = m_tuple9->addItem( "comp_costheta", m_comp_costheta );
    }
    else
    {
      log << MSG::ERROR << "    Cannot book N-tuple:" << long( m_tuple9 ) << endmsg;
      return StatusCode::FAILURE;
    }
  }
  // check pid and kaltrack
 
  //
  //--------End of book--------
  //
 
  log << MSG::INFO << "successfully return from initialize()" << endmsg;
  return StatusCode::SUCCESS;
}
 
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
StatusCode Single::execute() {
  setFilterPassed( false );
  // std::cout << "execute()" << std::endl;
 
  MsgStream log( msgSvc(), name() );
  log << MSG::INFO << "in execute()" << endmsg;
 
  // runinfo
  SmartDataPtr<Event::EventHeader> eventHeader( eventSvc(), "/Event/EventHeader" );
  int                              runNo = eventHeader->runNumber();
  int                              event = eventHeader->eventNumber();
  log << MSG::DEBUG << "run, evtnum = " << runNo << " , " << event << endmsg;
 
  // cout<<"event "<<event<<endl;
  Ncut0++;
 
  m_runNo = runNo;
  m_event = event;
 
  SmartDataPtr<EvtRecEvent> evtRecEvent( eventSvc(), EventModel::EvtRec::EvtRecEvent );
  //  log << MSG::INFO << "get event tag OK" << endmsg;
  log << MSG::DEBUG << "ncharg, nneu, tottks = " << evtRecEvent->totalCharged() << " , "
      << evtRecEvent->totalNeutral() << " , " << evtRecEvent->totalTracks() << endmsg;
 
  SmartDataPtr<EvtRecTrackCol> evtRecTrkCol( eventSvc(), EventModel::EvtRec::EvtRecTrackCol );
 
  Vint iGood;
  iGood.clear();
 
  int nCharge = 0;
 
  Hep3Vector    xorigin( 0, 0, 0 );
  IVertexDbSvc* vtxsvc;
  Gaudi::svcLocator()->service( "VertexDbSvc", vtxsvc );
  if ( vtxsvc->isVertexValid() )
  {
    double* dbv = vtxsvc->PrimaryVertex();
    double* vv  = vtxsvc->SigmaPrimaryVertex();
 
    xorigin.setX( dbv[0] );
    xorigin.setY( dbv[1] );
    xorigin.setZ( dbv[2] );
  }
 
  //
  // check McParticle infomation
  //
 
  SmartDataPtr<Event::McParticleCol> mcParticleCol( eventSvc(), "/Event/MC/McParticleCol" );
  if ( !mcParticleCol )
  {
    std::cout << "could not retrieve McParticleCol" << std::endl;
    return StatusCode::FAILURE;
  }
 
  double                         costhetamc = -9.;
  Event::McParticleCol::iterator iter_mc    = mcParticleCol->begin();
  HepLorentzVector               initialFourMomentum;
  for ( ; iter_mc != mcParticleCol->end(); iter_mc++ )
  {
    if ( abs( ( *iter_mc )->particleProperty() ) != abs( m_particleID ) ) continue;
    HepLorentzVector initialPosition = ( *iter_mc )->initialPosition();
    m_vx                             = initialPosition.x();
    m_vy                             = initialPosition.y();
    m_vz                             = initialPosition.z();
    m_vr0                            = sqrt( m_vx * m_vx + m_vy * m_vy );
 
    initialFourMomentum = ( *iter_mc )->initialFourMomentum();
    m_pxmc              = initialFourMomentum.px();
    m_pymc              = initialFourMomentum.py();
    m_pzmc              = initialFourMomentum.pz();
    m_e                 = initialFourMomentum.e();
    m_p                 = sqrt( m_pxmc * m_pxmc + m_pymc * m_pymc + m_pzmc * m_pzmc );
    m_pxymc             = sqrt( m_pxmc * m_pxmc + m_pymc * m_pymc );
    m_thetamc           = initialFourMomentum.theta() * 180 / ( CLHEP::pi );
    m_phimc             = initialFourMomentum.phi() * 180 / ( CLHEP::pi );
    m_m                 = initialFourMomentum.m();
    m_costhetamc        = initialFourMomentum.cosTheta();
 
    costhetamc = m_costhetamc;
 
    m_pro   = ( *iter_mc )->particleProperty();
    m_index = ( *iter_mc )->trackIndex();
  }
  m_mctotcharged = evtRecEvent->totalCharged();
 
  int mm = 0;
  if ( m_costhetamc > -0.7 && m_costhetamc < 0.7 && m_pxymc > 0.2 && m_pxymc < 0.3 ) mm = 1;
 
  double pxymc = m_pxymc;
  double pmc   = m_p;
  // if(m_mctotcharged==0&&m_costhetamc>-0.7&&m_costhetamc<-0.6&&m_pxymc>0.25&&m_pxymc<0.3)
  // cout<<"runNo="<<runNo<<"  event"<<event<<endl;
  //  if(m_mctotcharged==0&&m_costhetamc>-0.7&&m_costhetamc<-0.6&&m_pxymc>0.35&&m_pxymc<0.4)
  //  cout<<"runNo1="<<runNo<<"  event1"<<event<<endl;
  m_tuple8->write();
 
  //
  // check MdcTrack infomation
  //
 
  int itrack = 0;
 
  double Rz_comp   = 30.;
  double Rxy_comp  = 10.;
  double Rang_comp = 50.;
 
  for ( int i = 0; i < evtRecEvent->totalCharged(); i++ )
  {
    EvtRecTrackIterator itTrk = evtRecTrkCol->begin() + i;
    if ( !( *itTrk )->isMdcTrackValid() ) continue;
    if ( !( *itTrk )->isMdcKalTrackValid() ) continue;
    RecMdcTrack* mdcTrk    = ( *itTrk )->mdcTrack();
    RecMdcTrack* mdcKalTrk = ( *itTrk )->mdcTrack();
    double       pch       = mdcTrk->p();
    double       x0        = mdcTrk->x();
    double       y0        = mdcTrk->y();
    double       z0        = mdcTrk->z();
    int          trackid   = mdcTrk->trackId();
    int          charge    = mdcTrk->charge();
    double       pxy0      = mdcTrk->pxy();
    double       px0       = mdcTrk->px();
    double       py0       = mdcTrk->py();
    double       pz0       = mdcTrk->pz();
    double       r0        = mdcTrk->r();
    double       phi       = mdcTrk->phi() * 180 / ( CLHEP::pi );
    double       theta     = mdcTrk->theta() * 180 / ( CLHEP::pi );
    int          ndof      = mdcTrk->ndof();
    int          nster     = mdcTrk->nster();
    int          stat      = mdcTrk->stat();
 
    double phi0 = mdcTrk->helix( 1 );
    double xv   = xorigin.x();
    double yv   = xorigin.y();
    double Rxy  = ( x0 - xv ) * cos( phi0 ) + ( y0 - yv ) * sin( phi0 );
 
    m_vx0     = x0;
    m_vy0     = y0;
    m_vz0     = z0;
    m_pch     = pch;
    m_trackid = trackid;
    m_charge  = charge;
    m_pxy0    = pxy0;
    m_px0     = px0;
    m_py0     = py0;
    m_pz0     = pz0;
    m_r0      = r0;
    m_phi     = phi;
    m_theta   = theta;
    m_ndof    = ndof;
    m_nster   = nster;
    m_stat    = stat;
    m_vr0     = Rxy;
 
    HepVector a = mdcTrk->helix();
    m_d0        = a[0];
    m_phi0      = a[1];
    m_kappa     = a[2];
    m_dzhelix   = a[3];
    m_tanlamda  = a[4];
 
    HepSymMatrix Ea = mdcTrk->err();
 
    m_err11 = Ea[0][0];
    m_err21 = Ea[1][0];
    m_err22 = Ea[1][1];
    m_err31 = Ea[2][0];
    m_err32 = Ea[2][1];
    m_err33 = Ea[2][2];
    m_err41 = Ea[3][0];
    m_err42 = Ea[3][1];
    m_err43 = Ea[3][2];
    m_err44 = Ea[3][3];
    m_err51 = Ea[4][0];
    m_err52 = Ea[4][1];
    m_err53 = Ea[4][2];
    m_err54 = Ea[4][3];
    m_err55 = Ea[4][4];
 
    HepPoint3D point0( 0., 0., 0. ); // the initial point for MDC recosntruction
    HepPoint3D IP( xorigin[0], xorigin[1], xorigin[2] );
    VFHelix    helixip( point0, a, Ea );
    helixip.pivot( IP );
    HepVector vecipa = helixip.a();
    double    Rvxy0  = fabs( vecipa[0] ); // the nearest distance to IP in xy plane
    double    Rvz0   = vecipa[3];         // the nearest distance to IP in z direction
    double    Rvphi0 = vecipa[1];
    m_rvxy0          = Rvxy0;
    m_rvz0           = Rvz0;
    // m_rvphi0=Rvphi0;
 
    HepLorentzVector mdcpos( m_px0, m_py0, m_pz0, m_e );
    m_angd = initialFourMomentum.angle( mdcpos ) * 180 / ( CLHEP::pi );
 
    if ( fabs( Rvz0 ) >= m_vz0cut ) continue;
    if ( fabs( Rvxy0 ) >= m_vxy0cut ) continue;
    // if(m_angd>m_angdcut) continue;
 
    if ( m_angd < Rang_comp )
    {
      // if(Rvz0<Rz_comp||Rvxy0<Rxy_comp||m_angd<Rang_comp){
      Rz_comp   = Rvz0;
      Rxy_comp  = Rvxy0;
      Rang_comp = m_angd;
      itrack    = i;
    }
    //   nCharge += mdcTrk->charge();
    // if(nCharge!=-1) continue;
    // m_tuple1->write();
    iGood.push_back( i );
    nCharge += mdcTrk->charge();
  }
 
  // if(nCharge!=-1) return StatusCode::SUCCESS;
 
  //
  // Finish MdcTrack and Good Charged Track Selection
  //
 
  int nGood = iGood.size();
  log << MSG::DEBUG << "nGood, totcharge = " << nGood << " , " << nCharge << endmsg;
  // cout << "nGood, totcharge = " << nGood << " , " << nCharge << endl;
  if ( mm == 1 && nGood == 0 ) setFilterPassed( true );
 
  m_ngch    = nGood;
  m_ncharge = nCharge;
  for ( int i = 0; i < nGood; i++ )
  {
    EvtRecTrackIterator itTrk = evtRecTrkCol->begin() + iGood[i];
 
    RecMdcKalTrack* mdcKalTrk = ( *itTrk )->mdcKalTrack();
    if ( abs( m_particleID ) == 211 ) RecMdcKalTrack::setPidType( RecMdcKalTrack::pion );
    if ( abs( m_particleID ) == 321 ) RecMdcKalTrack::setPidType( RecMdcKalTrack::kaon );
    if ( abs( m_particleID ) == 2212 ) RecMdcKalTrack::setPidType( RecMdcKalTrack::proton );
    if ( abs( m_particleID ) == 11 ) RecMdcKalTrack::setPidType( RecMdcKalTrack::electron );
    if ( abs( m_particleID ) == 13 ) RecMdcKalTrack::setPidType( RecMdcKalTrack::muon );
 
    m_mdc_kal_px  = mdcKalTrk->px();
    m_mdc_kal_py  = mdcKalTrk->py();
    m_mdc_kal_pz  = mdcKalTrk->pz();
    m_mdc_kal_p   = sqrt( m_mdc_kal_px * m_mdc_kal_px + m_mdc_kal_py * m_mdc_kal_py +
                          m_mdc_kal_pz * m_mdc_kal_pz );
    m_mdc_kal_pxy = sqrt( m_mdc_kal_px * m_mdc_kal_px + m_mdc_kal_py * m_mdc_kal_py );
 
    m_mdc_kal_costheta = cos( mdcKalTrk->theta() );
  }
 
  double mdc_pxy0     = m_pxy0;
  double mdc_pch      = m_pch;
  double mdc_costheta = cos( m_theta );
 
  m_tuple1->write();
  // if(nGood!=1) return StatusCode::SUCCESS;
  // if(fabs(costhetamc)>0.93) return StatusCode::SUCCESS;
 
  Ncut1++;
 
  m_totcharged = evtRecEvent->totalCharged();
  m_totneutral = evtRecEvent->totalNeutral();
  m_tottracks  = evtRecEvent->totalTracks();
  m_ngood      = nGood;
 
  //  if(m_totneutral>0) return StatusCode::SUCCESS;
 
  // m_kaln=-9;
  // m_costhetan=-9;
 
  // for(int i = 0; i < nGood; i++) {
  // if (i != itrack) continue;
  m_kaln      = pxymc;
  m_costhetan = costhetamc;
  m_pn        = pmc;
  //}
  m_mat_charge    = -2;
  m_diff_pxy      = -2;
  m_diff_pch      = -2;
  m_diff_costheta = -2;
 
  for ( int i = 0; i < nGood; i++ )
  {
    if ( i != itrack ) continue;
    EvtRecTrackIterator itTrk  = evtRecTrkCol->begin() + iGood[i];
    RecMdcTrack*        mdcTrk = ( *itTrk )->mdcTrack();
    m_mat_charge               = mdcTrk->charge();
    m_diff_pxy                 = fabs( m_kaln - mdc_pxy0 );
    m_diff_pch                 = fabs( m_pn - mdc_pch );
    m_diff_costheta            = fabs( m_costhetan - mdc_costheta );
  }
 
  m_tuple7->write();
 
  // Finish runinfo check
 
  //
  //
  // check dedx infomation
  //
  //
  if ( m_checkDedx == 1 )
  {
    for ( int i = 0; i < nGood; i++ )
    {
      if ( m_itrack == 1 && i != itrack ) continue;
      EvtRecTrackIterator itTrk = evtRecTrkCol->begin() + i;
      if ( !( *itTrk )->isMdcTrackValid() ) continue;
      if ( !( *itTrk )->isMdcDedxValid() ) continue;
      RecMdcTrack* mdcTrk  = ( *itTrk )->mdcTrack();
      RecMdcDedx*  dedxTrk = ( *itTrk )->mdcDedx();
      m_ptrk               = mdcTrk->p();
 
      m_chie   = dedxTrk->chiE();
      m_chimu  = dedxTrk->chiMu();
      m_chipi  = dedxTrk->chiPi();
      m_chik   = dedxTrk->chiK();
      m_chip   = dedxTrk->chiP();
      m_ghit   = dedxTrk->numGoodHits();
      m_thit   = dedxTrk->numTotalHits();
      m_probPH = dedxTrk->probPH();
      m_normPH = dedxTrk->normPH();
      m_tuple2->write();
    }
  }
 
  //
  //
  // check ExtTrack  infomation
  //
  //
  if ( m_checkExt == 1 )
  {
    for ( int i = 0; i < nGood; i++ )
    {
      if ( m_itrack == 1 && i != itrack ) continue;
      EvtRecTrackIterator itTrk = evtRecTrkCol->begin() + i;
      if ( !( *itTrk )->isExtTrackValid() ) continue;
      RecExtTrack* extTrk = ( *itTrk )->extTrack();
 
      m_tof1     = extTrk->tof1();
      m_tof1path = extTrk->tof1Path();
      m_tof1z    = extTrk->tof1PosSigmaAlongZ();
      m_tof1t    = extTrk->tof1PosSigmaAlongT();
      m_tof1x    = extTrk->tof1PosSigmaAlongX();
      m_tof1y    = extTrk->tof1PosSigmaAlongY();
 
      m_tof2     = extTrk->tof2();
      m_tof2path = extTrk->tof2Path();
      m_tof2z    = extTrk->tof2PosSigmaAlongZ();
      m_tof2t    = extTrk->tof2PosSigmaAlongT();
      m_tof2x    = extTrk->tof2PosSigmaAlongX();
      m_tof2y    = extTrk->tof2PosSigmaAlongY();
 
      m_emctheta = extTrk->emcPosSigmaAlongTheta() * 180 / ( CLHEP::pi );
      m_emcphi   = extTrk->emcPosSigmaAlongPhi() * 180 / ( CLHEP::pi );
      m_emcpath  = extTrk->emcPath();
 
      m_mucz = extTrk->mucPosSigmaAlongZ();
      m_muct = extTrk->mucPosSigmaAlongT();
      m_mucx = extTrk->mucPosSigmaAlongX();
      m_mucy = extTrk->mucPosSigmaAlongY();
 
      m_tuple3->write();
    }
  }
 
  //
  //
  // check TOF  infomation
  //
  //
  if ( m_checkTof == 1 )
  {
    for ( int i = 0; i < nGood; i++ )
    {
      if ( m_itrack == 1 && i != itrack ) continue;
      EvtRecTrackIterator itTrk = evtRecTrkCol->begin() + i;
      if ( !( *itTrk )->isTofTrackValid() ) continue;
 
      SmartRefVector<RecTofTrack>           tofTrkCol = ( *itTrk )->tofTrack();
      SmartRefVector<RecTofTrack>::iterator iter_tof  = tofTrkCol.begin();
 
      m_path        = ( *iter_tof )->path();
      m_zrhit       = ( *iter_tof )->zrhit();
      m_ph          = ( *iter_tof )->ph();
      m_tof         = ( *iter_tof )->tof();
      m_errtof      = ( *iter_tof )->errtof();
      m_beta        = ( *iter_tof )->beta();
      m_texpe       = ( *iter_tof )->texpElectron();
      m_texpmu      = ( *iter_tof )->texpMuon();
      m_texppi      = ( *iter_tof )->texpPion();
      m_texpka      = ( *iter_tof )->texpKaon();
      m_texppro     = ( *iter_tof )->texpProton();
      m_toffsete    = ( *iter_tof )->toffsetElectron();
      m_toffsetmu   = ( *iter_tof )->toffsetMuon();
      m_toffsetpi   = ( *iter_tof )->toffsetPion();
      m_toffsetka   = ( *iter_tof )->toffsetKaon();
      m_toffsetpro  = ( *iter_tof )->toffsetProton();
      m_toffsetatpr = ( *iter_tof )->toffsetAntiProton();
      m_sigmae      = ( *iter_tof )->sigmaElectron();
      m_sigmamu     = ( *iter_tof )->sigmaMuon();
      m_sigmapi     = ( *iter_tof )->sigmaPion();
      m_sigmaka     = ( *iter_tof )->sigmaKaon();
      m_sigmapro    = ( *iter_tof )->sigmaProton();
      m_sigmaatpr   = ( *iter_tof )->sigmaAntiProton();
      m_quality     = ( *iter_tof )->quality();
      m_t0          = ( *iter_tof )->t0();
      m_errt0       = ( *iter_tof )->errt0();
      m_errz        = ( *iter_tof )->errz();
      m_phitof      = ( *iter_tof )->phi() * 180 / ( CLHEP::pi );
      m_errphi      = ( *iter_tof )->errphi() * 180 / ( CLHEP::pi );
      m_energy      = ( *iter_tof )->energy();
      m_errenergy   = ( *iter_tof )->errenergy();
 
      m_tuple4->write();
    }
  }
 
  //
  //
  // check EmcTrack  infomation
  //
  //
  if ( m_checkEmc == 1 )
  {
    for ( int i = 0; i < nGood; i++ )
    {
      if ( m_itrack == 1 && i != itrack ) continue;
      EvtRecTrackIterator itTrk = evtRecTrkCol->begin() + i;
      if ( !( *itTrk )->isEmcShowerValid() ) continue;
      RecEmcShower* emcTrk = ( *itTrk )->emcShower();
 
      m_x         = emcTrk->x();
      m_y         = emcTrk->y();
      m_z         = emcTrk->z();
      m_thetaemc  = emcTrk->theta() * 180 / ( CLHEP::pi );
      m_phiemc    = emcTrk->phi() * 180 / ( CLHEP::pi );
      m_dx        = emcTrk->dx();
      m_dy        = emcTrk->dy();
      m_dz        = emcTrk->dz();
      m_dtheta    = emcTrk->dtheta() * 180 / ( CLHEP::pi );
      m_dphi      = emcTrk->dphi() * 180 / ( CLHEP::pi );
      m_energyemc = emcTrk->energy();
      m_de        = emcTrk->dE();
      m_eseed     = emcTrk->eSeed();
      m_e3x3      = emcTrk->e3x3();
      m_e5x5      = emcTrk->e5x5();
      m_secp      = emcTrk->secondMoment();
      m_latp      = emcTrk->latMoment();
 
      m_tuple5->write();
    }
  }
 
  //
  // Assign 4-momentum to each charged track
  //
  ParticleID* pid = ParticleID::instance();
  for ( int i = 0; i < nGood; i++ )
  {
    if ( nGood != 1 ) continue;
    if ( m_itrack == 1 && i != itrack ) continue;
    EvtRecTrackIterator itTrk = evtRecTrkCol->begin() + iGood[i];
    //    if(pid) delete pid;
 
    pid->init();
    pid->setMethod( pid->methodProbability() );
    // pid->setMethod(pid->methodLikelihood());  //for Likelihood Method
 
    pid->setChiMinCut( 4 );
    pid->setRecTrack( *itTrk );
    if ( m_pid_mode == 0 )
      pid->usePidSys( pid->useDedx() | pid->useTof1() | pid->useTof2() | pid->useTofE() );
    if ( m_pid_mode == 1 ) pid->usePidSys( pid->useDedx() );
    if ( m_pid_mode == 2 ) pid->usePidSys( pid->useTof1() | pid->useTof2() | pid->useTofE() );
 
    // pid->identify(pid->onlyPion() | pid->onlyKaon());    // seperater Pion/Kaon
    //   pid->identify(pid->onlyPion());
    pid->identify( pid->onlyPionKaonProton() );
    pid->calculate();
 
    // if(!(pid->IsPidInfoValid())&&costhetamc>-0.7&&costhetamc<-0.6&&pxymc>0.25&&pxymc<0.3)
    // cout<<"pidrunNo="<<runNo<<"  pidevent"<<event<<endl;
    // if(!(pid->IsPidInfoValid())&&costhetamc>-0.7&&costhetamc<-0.6&&pxymc>0.35&&pxymc<0.4)
    // cout<<"pidrunNo1="<<runNo<<" pidevent1"<<event<<endl;
    if ( !( pid->IsPidInfoValid() ) ) continue;
    RecMdcTrack* mdcTrk = ( *itTrk )->mdcTrack();
    m_ptrk_pid          = mdcTrk->p();
    m_cost_pid          = cos( mdcTrk->theta() );
    m_dedx_pid          = pid->chiDedx( 2 );
    m_tof1_pid          = pid->chiTof1( 2 );
    m_tof2_pid          = pid->chiTof2( 2 );
    m_prob_pid_pion     = pid->probPion();
    m_prob_pid_kaon     = pid->probKaon();
    m_prob_pid_proton   = pid->probProton();
 
    //
    // check KalTrack  infomation
    //
 
    // if(!(*itTrk)->isMdcKalTrackValid()) {cout<<"kal"<<endl; continue;}
 
    RecMdcKalTrack* mdcKalTrk = ( *itTrk )->mdcKalTrack();
    if ( abs( m_particleID ) == 211 ) RecMdcKalTrack::setPidType( RecMdcKalTrack::pion );
    if ( abs( m_particleID ) == 321 ) RecMdcKalTrack::setPidType( RecMdcKalTrack::kaon );
    if ( abs( m_particleID ) == 2212 ) RecMdcKalTrack::setPidType( RecMdcKalTrack::proton );
    if ( abs( m_particleID ) == 11 ) RecMdcKalTrack::setPidType( RecMdcKalTrack::electron );
    if ( abs( m_particleID ) == 13 ) RecMdcKalTrack::setPidType( RecMdcKalTrack::muon );
 
    m_kalpx  = mdcKalTrk->px();
    m_kalpy  = mdcKalTrk->py();
    m_kalpz  = mdcKalTrk->pz();
    m_kalp   = sqrt( m_kalpx * m_kalpx + m_kalpy * m_kalpy + m_kalpz * m_kalpz );
    m_kalpxy = sqrt( m_kalpx * m_kalpx + m_kalpy * m_kalpy );
 
    HepLorentzVector ptr;
 
    ptr           = mdcKalTrk->p4( 0.000511 * 0.000511 );
    m_kalcostheta = ptr.cosTheta();
 
    double kalcostheta = m_kalcostheta;
    m_comp_costheta    = costhetamc - kalcostheta;
    // cout<<"   kalcostheta="<<m_kalcostheta<<endl;
    HepVector    k  = mdcKalTrk->helix();
    HepSymMatrix Ek = mdcKalTrk->err();
 
    m_d0kal       = k[0];
    m_phi0kal     = k[1];
    m_kappakal    = k[2];
    m_dzhelixkal  = k[3];
    m_tanlamdakal = k[4];
 
    m_err11kal = Ek[0][0];
    m_err21kal = Ek[1][0];
    m_err22kal = Ek[1][1];
 
    m_err31kal = Ek[2][0];
    m_err32kal = Ek[2][1];
    m_err33kal = Ek[2][2];
    m_err41kal = Ek[3][0];
    m_err42kal = Ek[3][1];
    m_err43kal = Ek[3][2];
    m_err44kal = Ek[3][3];
    m_err51kal = Ek[4][0];
    m_err52kal = Ek[4][1];
    m_err53kal = Ek[4][2];
    m_err54kal = Ek[4][3];
    m_err55kal = Ek[4][4];
 
    m_kalnn      = pxymc;
    m_costhetann = costhetamc;
    m_pnn        = pmc;
    // cout<<"pxy="<<pxymc<<'\t'<<'\t'<<m_kalpxy<<endl;
    // cout<<"costheta="<<costhetamc<<'\t'<<m_kalcostheta<<endl<<endl;
 
    m_tuple9->write();
  }
 
  //
  //
  // check MucTrack  infomation
  //
  //
  if ( m_checkMuc == 1 )
  {
    for ( int i = 0; i < nGood; i++ )
    {
      if ( m_itrack == 1 && i != itrack ) continue;
      EvtRecTrackIterator itTrk = evtRecTrkCol->begin() + i;
      if ( !( *itTrk )->isMucTrackValid() ) continue;
      RecMucTrack* mucTrk = ( *itTrk )->mucTrack();
 
      m_depth     = mucTrk->depth();
      m_chi2      = mucTrk->chi2();
      m_rms       = mucTrk->rms();
      m_xpos      = mucTrk->xPos();
      m_ypos      = mucTrk->yPos();
      m_zpos      = mucTrk->zPos();
      m_xpossigma = mucTrk->xPosSigma();
      m_ypossigma = mucTrk->yPosSigma();
      m_zpossigma = mucTrk->zPosSigma();
      m_px        = mucTrk->px();
      m_py        = mucTrk->py();
      m_pz        = mucTrk->pz();
      m_distance  = mucTrk->distance();
      m_deltaphi  = mucTrk->deltaPhi() * 180 / ( CLHEP::pi );
 
      m_tuple6->write();
    }
  }
 
  return StatusCode::SUCCESS;
}
 
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
StatusCode Single::finalize() {
  cout << "mdc" << endl;
  cout << "total number:         " << Ncut0 << endl;
  cout << "nGood: " << Ncut1 << endl;
  MsgStream log( msgSvc(), name() );
  log << MSG::INFO << "in finalize()" << endmsg;
  return StatusCode::SUCCESS;
}