BesMucEfficiency.cc

Go to the documentation of this file.

//
//
//
//
 
#include "GaudiKernel/Bootstrap.h"
#include "GaudiKernel/ISvcLocator.h"
#include "GaudiKernel/MsgStream.h"
 
#include "G4Box.hh"
#include "G4VSolid.hh"
#include "TFile.h"
#include "TTree.h"
#include <fstream>
#include <iostream>
#include <sstream>
 
#include "MucSim/BesMucEfficiency.hh"
 
using namespace std;
BesMucEfficiency* BesMucEfficiency::fPointer = 0;
BesMucEfficiency* BesMucEfficiency::Instance( void ) {
  if ( !fPointer ) fPointer = new BesMucEfficiency();
  return fPointer;
}
 
BesMucEfficiency::BesMucEfficiency() {
  if ( fPointer )
  {
    G4cout << "BesMucEfficiency constructed twice." << G4endl;
    exit( -1 );
  }
  fPointer = this;
}
 
BesMucEfficiency::~BesMucEfficiency() {}
 
void BesMucEfficiency::Initialize( G4String filename ) {
  G4int    part, seg, gap, strip, pad;
  G4double effi;
  //  for(G4int part=0;part<3;part++){
  //    for(G4int seg=0;seg<8;seg++){
  //      for(G4int gap=0;gap<gap_Max;gap++){
  //    for(G4int strip=0;strip<strip_Max;strip++){
  //      for(G4int pad=0;pad<pad_Max;pad++){
  //            m_effi[part][seg][gap][strip][pad]=1;
  //      }
  //    }
  //    }
  //  }
  //}
  // G4cout<<"in BesMucEfficiency::Initialize()"<<G4endl;
 
  // TFile *f=new TFile("muc-effi.root");
  // TTree *t1=(TTree*)f->Get("t1");
 
  std::ifstream fin( filename );
 
  char  buffer[100];
  G4int num = 0;
 
  fin.getline( buffer, 100, '\n' ); // get info whether add effi or not
  std::istringstream stringBuf( buffer );
  stringBuf >> IsAddEffi;
  // G4cout<<"IsAddEffi ="<<IsAddEffi<<G4endl;
  fin.getline( buffer, 100, '\n' );
 
  if ( !fin )
  {
    G4cout << "error opening effi data" << G4endl;
    IsAddEffi = 1.0; // no effi data. set effi = 1.0
  }
 
  //      fin.getline(buffer,100,'\n');
  //      std::istringstream stringBuf2(buffer);
  //      stringBuf2>>part>>seg>>gap>>strip>>pad;
  //      G4cout<<"---------- "<<pad<<endl;
 
  while ( fin.getline( buffer, 100, '\n' ) )
  {
    std::istringstream stringBuf2( buffer );
    stringBuf2 >> part >> seg >> gap >> strip >> pad >> effi;
    m_effi[part][seg][gap][strip][pad] = effi;
    num++;
  }
  for ( G4int seg = 0; seg < 8; seg++ )
  {
    for ( G4int strip = 0; strip < strip_Max; strip++ ) { m_effi[1][seg][0][strip][105] = 0; }
  }
 
  // G4cout<<"------------in Effi::init()----  "<<num<<G4endl;
}
 
void BesMucEfficiency::CheckCalibSvc() {
 
  ISvcLocator* svcLocator = Gaudi::svcLocator();
  // IMucCalibConstSvc* m_ptrCalibSvc;
  StatusCode sc = svcLocator->service( "MucCalibConstSvc", m_ptrCalibSvc, true );
 
  if ( sc != StatusCode::SUCCESS ) { G4cout << "Can not use MucCalibConstSvc!" << G4endl; }
}
 
void BesMucEfficiency::SetHit( BesMucHit* hit ) {
  m_pHit     = hit;
  G4int part = m_pHit->GetPart();
  G4int gap  = m_pHit->GetGap();
  m_Pos_Hit  = m_pHit->GetPosLocal().y(); // different from BesMucdigit
  if ( ( part == 1 && gap % 2 != 0 ) || ( part != 1 && gap % 2 == 0 ) )
  { m_Pos_Hit = m_pHit->GetPosLocal().x(); }
 
  // set m_Strip, m_Pos_Strip, m_Length, m_Width
  BesMucDigit aDigit;
  aDigit.SetHit( m_pHit );
  m_Strip = aDigit.GetNearestStripNo();
  //  G4VPhysicalVolume* pvGasChamber = m_pHit->GetVolume();
  GetPosLengthWidth( aDigit.GetNearestStrip() );
 
  //  G4cout<<"m_Pos_Hit = "<<m_Pos_Hit<<G4endl;
}
 
G4int BesMucEfficiency::GetPad() { // it will be better to put this function into Class
                                   // BesMucDigit
  G4double pad1 = ( m_Pos_Hit + m_Length / 2 - m_Pos_Strip ) / m_Width;
  G4int    pad  = G4int( pad1 );
  // G4cout<<"---in Effi::GetPad()--- hit: "<<m_Pos_Hit<<" part: "<<m_pHit->GetPart()<<" gap:
  // "<<m_pHit->GetGap()<<"  L: "<<m_Length<<"  strip: "<<m_Pos_Strip<<"  width: "<<m_Width<<"
  // pad: "<<pad<<G4endl;
  if ( abs( m_Pos_Hit - m_Pos_Strip ) < m_Length / 2 ) return pad;
  else return -999;
}
 
G4double BesMucEfficiency::GetEfficiency() {
  // look up table with (part;seg;gap;m_Strip)
  G4int part  = m_pHit->GetPart();
  G4int seg   = m_pHit->GetSeg();
  G4int gap   = m_pHit->GetGap();
  G4int strip = m_Strip;
  G4int pad   = GetPad();
 
  G4double eff = 0;
  if ( 0 != m_ptrCalibSvc )
  {
    eff = m_ptrCalibSvc->getEff( part, seg, gap, strip );
    // G4cout << "Prt: " << part << "\tseg: " << seg << "\tlay: " << gap << "\tstr: " <<
    // m_Strip
    //           << "\t:" << eff << endl;
  }
  else
  {
    // G4cout << "CalibSvc unavailable!" << G4endl;
    eff = 0.95;
  }
  // G4cout<<part<<"\t"<<seg<<"\t"<<gap<<"\t"<<m_Strip<<"\t"<<eff<<G4endl;
  return eff;
}
 
void BesMucEfficiency::GetPosLengthWidth( G4VPhysicalVolume* pvStrip ) {
  m_Pos_Strip = 1.0e38;
 
  G4int part = m_pHit->GetPart();
  G4int gap  = m_pHit->GetGap();
 
  m_Pos_Strip = pvStrip->GetObjectTranslation().y();
  if ( ( part == 1 && gap % 2 != 0 ) || ( part != 1 && gap % 2 == 0 ) )
  { m_Pos_Strip = pvStrip->GetObjectTranslation().x(); }
 
  G4String striptype = pvStrip->GetLogicalVolume()->GetSolid()->GetEntityType();
  //  G4String striplenght= pvStrip->GetLogicalVolume()->GetName();
  G4Box* temp;
  temp     = (G4Box*)pvStrip->GetLogicalVolume()->GetSolid();
  m_Width  = temp->GetXHalfLength() * 2;
  m_Length = temp->GetYHalfLength() * 2;
  if ( ( part == 1 && gap % 2 != 0 ) || ( part != 1 && gap % 2 == 0 ) )
  {
    m_Width  = temp->GetYHalfLength() * 2;
    m_Length = temp->GetXHalfLength() * 2;
  }
  // G4cout<<"in Set  "<<m_Length<<" "<<temp->GetXHalfLength()<<" "<<m_Width<<"
  // "<<m_Pos_Strip<<G4endl;
}