BOSS 6.6.4.p01
BESIII Offline Software System
Loading...
Searching...
No Matches
BesMucEfficiency Class Reference

#include <BesMucEfficiency.hh>

Public Member Functions

 BesMucEfficiency ()
 
 ~BesMucEfficiency ()
 
void Initialize (G4String filename)
 
void CheckCalibSvc ()
 
void SetHit (BesMucHit *hit)
 
void GetPosLengthWidth (G4VPhysicalVolume *pvStrip)
 
G4int GetPad ()
 
G4double GetEfficiency ()
 

Static Public Member Functions

static BesMucEfficiencyInstance (void)
 

Public Attributes

IMessageSvc * msgSvc
 

Detailed Description

Definition at line 31 of file BesMucEfficiency.hh.

Constructor & Destructor Documentation

◆ BesMucEfficiency()

BesMucEfficiency::BesMucEfficiency ( )

Definition at line 27 of file BesMucEfficiency.cc.

28{
29 if(fPointer)
30 {G4Exception("BesMucEfficiency constructed twice.");}
31 fPointer=this;
32
33}

Referenced by Instance().

◆ ~BesMucEfficiency()

BesMucEfficiency::~BesMucEfficiency ( )

Definition at line 35 of file BesMucEfficiency.cc.

36{
37}

Member Function Documentation

◆ CheckCalibSvc()

void BesMucEfficiency::CheckCalibSvc ( )

Definition at line 96 of file BesMucEfficiency.cc.

97{
98
99 ISvcLocator* svcLocator = Gaudi::svcLocator();
100 //IMucCalibConstSvc* m_ptrCalibSvc;
101 StatusCode sc = svcLocator->service("MucCalibConstSvc", m_ptrCalibSvc, true);
102
103 if( sc != StatusCode::SUCCESS){
104 G4cout<< "Can not use MucCalibConstSvc!" << G4endl;
105 }
106
107}

Referenced by BesMucSD::ProcessHits().

◆ GetEfficiency()

G4double BesMucEfficiency::GetEfficiency ( )

Definition at line 141 of file BesMucEfficiency.cc.

142{
143 // look up table with (part;seg;gap;m_Strip)
144 G4int part = m_pHit->GetPart();
145 G4int seg = m_pHit->GetSeg();
146 G4int gap = m_pHit->GetGap();
147 G4int strip = m_Strip;
148 G4int pad = GetPad();
149
150 G4double eff = 0;
151 if( 0 != m_ptrCalibSvc ){
152 eff = m_ptrCalibSvc->getEff(part, seg, gap, strip);
153 //G4cout << "Prt: " << part << "\tseg: " << seg << "\tlay: " << gap << "\tstr: " << m_Strip
154 // << "\t:" << eff << endl;
155 }
156 else
157 {
158 //G4cout << "CalibSvc unavailable!" << G4endl;
159 eff = 0.95;
160 }
161 //G4cout<<part<<"\t"<<seg<<"\t"<<gap<<"\t"<<m_Strip<<"\t"<<eff<<G4endl;
162 return eff;
163}
G4int GetSeg()
Definition: BesMucHit.hh:75
G4int GetPart()
Definition: BesMucHit.hh:74
G4int GetGap()
Definition: BesMucHit.hh:76
virtual double getEff(int part, int segment, int layer, int strip) const =0

Referenced by BesMucSD::ProcessHits().

◆ GetPad()

G4int BesMucEfficiency::GetPad ( )

Definition at line 129 of file BesMucEfficiency.cc.

130{//it will be better to put this function into Class BesMucDigit
131 G4double pad1 = (m_Pos_Hit+m_Length/2-m_Pos_Strip)/m_Width;
132 G4int pad =G4int(pad1);
133 //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;
134 if(abs(m_Pos_Hit-m_Pos_Strip)<m_Length/2)
135 return pad;
136 else
137 return -999;
138
139}
double abs(const EvtComplex &c)
Definition: EvtComplex.hh:212

Referenced by GetEfficiency().

◆ GetPosLengthWidth()

void BesMucEfficiency::GetPosLengthWidth ( G4VPhysicalVolume *  pvStrip)

Definition at line 165 of file BesMucEfficiency.cc.

166{
167 m_Pos_Strip = 1.0e38;
168
169 G4int part = m_pHit->GetPart();
170 G4int gap = m_pHit->GetGap();
171
172 m_Pos_Strip = pvStrip->GetObjectTranslation().y();
173 if ( (part == 1 && gap%2 != 0) || (part != 1 && gap%2 == 0) ) {
174 m_Pos_Strip = pvStrip->GetObjectTranslation().x();
175 }
176
177 G4String striptype= pvStrip->GetLogicalVolume()->GetSolid()->GetEntityType();
178// G4String striplenght= pvStrip->GetLogicalVolume()->GetName();
179 G4Box *temp;
180 temp = (G4Box *)pvStrip->GetLogicalVolume()->GetSolid();
181 m_Width = temp->GetXHalfLength()*2;
182 m_Length = temp->GetYHalfLength()*2;
183 if ( (part == 1 && gap%2 != 0) || (part != 1 && gap%2 == 0) ) {
184 m_Width =temp->GetYHalfLength()*2;
185 m_Length=temp->GetXHalfLength()*2;
186 }
187 //G4cout<<"in Set "<<m_Length<<" "<<temp->GetXHalfLength()<<" "<<m_Width<<" "<<m_Pos_Strip<<G4endl;
188
189}

Referenced by SetHit().

◆ Initialize()

void BesMucEfficiency::Initialize ( G4String  filename)

Definition at line 40 of file BesMucEfficiency.cc.

41{ G4int part,seg,gap,strip,pad;
42 G4double effi;
43// for(G4int part=0;part<3;part++){
44// for(G4int seg=0;seg<8;seg++){
45// for(G4int gap=0;gap<gap_Max;gap++){
46// for(G4int strip=0;strip<strip_Max;strip++){
47// for(G4int pad=0;pad<pad_Max;pad++){
48// m_effi[part][seg][gap][strip][pad]=1;
49// }
50// }
51 // }
52 // }
53 //}
54 // G4cout<<"in BesMucEfficiency::Initialize()"<<G4endl;
55
56 // TFile *f=new TFile("muc-effi.root");
57 // TTree *t1=(TTree*)f->Get("t1");
58
59
60 std::ifstream fin(filename);
61
62 char buffer[100];
63 G4int num=0;
64
65 fin.getline(buffer,100,'\n'); //get info whether add effi or not
66 std::istringstream stringBuf(buffer);
67 stringBuf>>IsAddEffi;
68 //G4cout<<"IsAddEffi ="<<IsAddEffi<<G4endl;
69 fin.getline(buffer,100,'\n');
70
71 if(!fin){
72 G4cout<<"error opening effi data"<<G4endl;
73 IsAddEffi = 1.0; // no effi data. set effi = 1.0
74 }
75
76// fin.getline(buffer,100,'\n');
77// std::istringstream stringBuf2(buffer);
78// stringBuf2>>part>>seg>>gap>>strip>>pad;
79// G4cout<<"---------- "<<pad<<endl;
80
81 while(fin.getline(buffer,100,'\n')){
82 std::istringstream stringBuf2(buffer);
83 stringBuf2>>part>>seg>>gap>>strip>>pad>>effi;
84 m_effi[part][seg][gap][strip][pad] = effi;
85 num++;
86 }
87 for(G4int seg=0;seg<8;seg++){
88 for(G4int strip=0;strip<strip_Max;strip++){
89 m_effi[1][seg][0][strip][105]=0;
90 }
91 }
92
93 //G4cout<<"------------in Effi::init()---- "<<num<<G4endl;
94}
#define strip_Max
int num[96]
Definition: ranlxd.c:373

◆ Instance()

BesMucEfficiency * BesMucEfficiency::Instance ( void  )
static

Definition at line 21 of file BesMucEfficiency.cc.

21 {
22 if(!fPointer)fPointer = new BesMucEfficiency();
23 return fPointer;
24
25}

Referenced by BesMucSD::BesMucSD().

◆ SetHit()

void BesMucEfficiency::SetHit ( BesMucHit hit)

Definition at line 109 of file BesMucEfficiency.cc.

110{
111 m_pHit = hit;
112 G4int part = m_pHit->GetPart();
113 G4int gap = m_pHit->GetGap();
114 m_Pos_Hit = m_pHit->GetPosLocal().y(); //different from BesMucdigit
115 if ( (part == 1 && gap%2 != 0) || (part != 1 && gap%2 == 0) ) {
116 m_Pos_Hit = m_pHit->GetPosLocal().x();
117 }
118
119 //set m_Strip, m_Pos_Strip, m_Length, m_Width
120 BesMucDigit aDigit;
121 aDigit.SetHit(m_pHit);
122 m_Strip = aDigit.GetNearestStripNo();
123 // G4VPhysicalVolume* pvGasChamber = m_pHit->GetVolume();
125
126// G4cout<<"m_Pos_Hit = "<<m_Pos_Hit<<G4endl;
127}
G4int GetNearestStripNo()
Definition: BesMucDigit.cc:63
G4VPhysicalVolume * GetNearestStrip()
Definition: BesMucDigit.cc:53
void SetHit(BesMucHit *hit)
Definition: BesMucDigit.cc:29
void GetPosLengthWidth(G4VPhysicalVolume *pvStrip)
G4ThreeVector GetPosLocal()
Definition: BesMucHit.hh:69

Referenced by BesMucSD::ProcessHits().

Member Data Documentation

◆ msgSvc

IMessageSvc* BesMucEfficiency::msgSvc

Definition at line 48 of file BesMucEfficiency.hh.


The documentation for this class was generated from the following files: