BOSS 8.0.0
BESIII Offline Software System
Loading...
Searching...
No Matches
calib_endcap_atten Class Reference

#include <calib_endcap_atten.h>

Inheritance diagram for calib_endcap_atten:

Public Member Functions

 calib_endcap_atten (const unsigned int nrbin)
 ~calib_endcap_atten ()
void calculate (RecordSet *&data, unsigned int icounter)
Public Member Functions inherited from TofCalibFit
 TofCalibFit (bool isbarrel, const int npar)
 ~TofCalibFit ()
const std::string & name () const
void fillTxt (const char *file)
void fillRoot (const char *file)
HepVector tcorrelation ()
void setTCorrelation (HepVector tc)

Additional Inherited Members

Protected Attributes inherited from TofCalibFit
int m_npar
unsigned int nKind
unsigned int nBinPerCounter
unsigned int nHistPerCounter
unsigned int nCanvasPerCounter
std::vector< unsigned int > nGraphPerCanvasPerCounter
unsigned int nHistogram
unsigned int nCanvas
std::vector< unsigned int > nGraphPerCanvas
string m_name
HepVector X
std::vector< TH1F * > m_histograms
std::vector< TH1F * > m_graphs
std::vector< HepVector > m_result
std::vector< string > CanvasPerCounterName
std::vector< string > CanvasName
HepVector m_tcorrelation

Detailed Description

Definition at line 10 of file calib_endcap_atten.h.

Constructor & Destructor Documentation

◆ calib_endcap_atten()

calib_endcap_atten::calib_endcap_atten ( const unsigned int nrbin)

Definition at line 12 of file calib_endcap_atten.cxx.

13 : TofCalibFit( false, nEndcapAtten ) {
14
15 nKind = 1; // pulse height
16 nBinPerCounter = nrbin + 1;
17
20 CanvasPerCounterName.push_back( static_cast<string>( "Pulse Height Most Probable Value" ) );
21 CanvasPerCounterName.push_back( static_cast<string>( "Pulse Height Sigma" ) );
22 nGraphPerCanvasPerCounter.push_back( 1 );
23 nGraphPerCanvasPerCounter.push_back( 1 );
24
25 nHistogram = 0;
26 nCanvas = 2;
27 CanvasName.push_back(
28 static_cast<string>( "Pulse Height Most Probable Value vs TOF Counter Number" ) );
29 CanvasName.push_back( static_cast<string>( "Pulse Height Sigma vs TOF Counter Number" ) );
30 nGraphPerCanvas.push_back( 1 );
31 nGraphPerCanvas.push_back( 1 );
32
33 int numGraphs1 = 0;
34 std::vector<unsigned int>::iterator iter = nGraphPerCanvasPerCounter.begin();
35 for ( ; iter != nGraphPerCanvasPerCounter.end(); iter++ )
36 { numGraphs1 = numGraphs1 + ( *iter ); }
37 if ( numGraphs1 != nGraphEcAtten )
38 {
39 cout << "tofcalgsec::calib_endcap_atten: the number of Graphs is NOT reasonable!!!"
40 << endl;
41 exit( 0 );
42 }
43 int numGraphs2 = 0;
44 iter = nGraphPerCanvas.begin();
45 for ( ; iter != nGraphPerCanvas.end(); iter++ ) { numGraphs2 = numGraphs2 + ( *iter ); }
46 if ( numGraphs2 != nGraphEcAtten )
47 {
48 cout << "tofcalgsec::calib_endcap_atten: the number of Graphs is NOT reasonable!!!"
49 << endl;
50 exit( 0 );
51 }
52
53 m_name = string( "calib_endcap_atten" );
54
55 const int qbin = 100;
56 const double qbegin = 0.0;
57 const double qend = 5000.0;
58
59 // histograms per counter
60 char hname[256];
61 for ( unsigned int i = 0; i < NEndcap; i++ )
62 {
63 m_result.push_back( HepVector( nEndcapAtten, 0 ) );
64 for ( unsigned int k = 0; k < nrbin; k++ )
65 {
66 sprintf( hname, "Q-id%i-r%i", i, k );
67 m_histograms.push_back( new TH1F( hname, hname, qbin, qbegin, qend ) );
68 m_fitresult.push_back( HepVector( nParEcAtten, 0 ) );
69 }
70 sprintf( hname, "Q0-id%i", i );
71 m_histograms.push_back( new TH1F( hname, hname, qbin, qbegin, qend ) );
72 m_fitresult.push_back( HepVector( nParEcAtten, 0 ) );
73 }
74
75 rpos.resize( nrbin );
76 rposerr.resize( nrbin );
77 rstep = ( rtofend - rtofbegin ) / nrbin;
78 for ( unsigned int i = 0; i < nrbin; i++ )
79 {
80 rpos[i] = rtofbegin + ( i + 0.5 ) * rstep;
81 rposerr[i] = 0.5 * rstep;
82 }
83 itofid.resize( NEndcap );
84 itofiderr.resize( NEndcap );
85 itofidstep = 1.0;
86 for ( unsigned int i = 0; i < NEndcap; i++ )
87 {
88 itofid[i] = i * 1.0;
89 itofiderr[i] = 0.5;
90 }
91}
sprintf(cut, "kal_costheta0_em>-0.93&&kal_costheta0_em<0.93&&kal_pxy0_em>=0.05+%d*0.1&&kal_" "pxy0_em<0.15+%d*0.1&&NGch>=2", j, j)
EvtStreamInputIterator< typename Generator::result_type > iter(Generator gen, int N=0)
const double rtofend
Definition TofCalibFit.h:17
const double rtofbegin
Definition TofCalibFit.h:16
const unsigned int NEndcap
Definition TofDataSet.h:13
const int nEndcapAtten
const int nGraphEcAtten
const int nParEcAtten
string m_name
Definition TofCalibFit.h:50
std::vector< HepVector > m_result
Definition TofCalibFit.h:55
unsigned int nCanvas
Definition TofCalibFit.h:47
std::vector< TH1F * > m_histograms
Definition TofCalibFit.h:53
unsigned int nBinPerCounter
Definition TofCalibFit.h:41
unsigned int nKind
Definition TofCalibFit.h:40
std::vector< unsigned int > nGraphPerCanvasPerCounter
Definition TofCalibFit.h:45
unsigned int nCanvasPerCounter
Definition TofCalibFit.h:44
unsigned int nHistPerCounter
Definition TofCalibFit.h:43
std::vector< unsigned int > nGraphPerCanvas
Definition TofCalibFit.h:48
std::vector< string > CanvasName
Definition TofCalibFit.h:58
unsigned int nHistogram
Definition TofCalibFit.h:46
TofCalibFit(bool isbarrel, const int npar)
std::vector< string > CanvasPerCounterName
Definition TofCalibFit.h:57

◆ ~calib_endcap_atten()

calib_endcap_atten::~calib_endcap_atten ( )

Definition at line 93 of file calib_endcap_atten.cxx.

93 {
94 m_fitresult.clear();
95 rpos.clear();
96 rposerr.clear();
97 itofid.clear();
98 itofiderr.clear();
99}

Member Function Documentation

◆ calculate()

void calib_endcap_atten::calculate ( RecordSet *& data,
unsigned int icounter )
virtual

Implements TofCalibFit.

Definition at line 101 of file calib_endcap_atten.cxx.

101 {
102
103 std::cout << setiosflags( ios::left ) << setw( 10 ) << icounter << setw( 8 ) << data->size()
104 << setw( 30 ) << name() << std::endl;
105
106 if ( data->size() > 0 )
107 {
108 std::vector<Record*>::iterator iter = data->begin();
109 for ( ; iter != data->end(); iter++ ) { fillRecord( ( *iter ), icounter ); }
110 fitHistogram( icounter );
111 fillGraph( icounter );
112 fitGraph( icounter );
113 }
114 else
115 {
116 fillGraph( icounter ); // keep the m_graphs is not empty()
117 }
118
119 if ( data->size() > 0 )
120 {
121 std::vector<Record*>::iterator iter = data->begin();
122 for ( ; iter != data->end(); iter++ )
123 {
124 updateData( ( *iter ), icounter );
125 fillRecordQ0( ( *iter ), icounter );
126 }
127 fitHistogramQ0( icounter );
128 }
129
130 if ( icounter == ( NEndcap - 1 ) ) { fillGraphQ0(); }
131
132 return;
133}
TTree * data
const std::string & name() const
Definition TofCalibFit.h:27

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