#include <PmtCalibLeadingEdgeWithCuts.h>
Inheritance diagram for PmtCalibLeadingEdgeWithCuts:
Public Types | |
| SUCCESS | |
| NO_INTERFACE | |
| VERSMISMATCH | |
| LAST_ERROR | |
| SUCCESS | |
| NO_INTERFACE | |
| VERSMISMATCH | |
| LAST_ERROR | |
| enum | Status |
Public Member Functions | |
| PmtCalibLeadingEdgeWithCuts (const std::string &type, const std::string &name, const IInterface *parent) | |
| virtual | ~PmtCalibLeadingEdgeWithCuts () |
| virtual StatusCode | initialize () |
| virtual StatusCode | finalize () |
| virtual StatusCode | process (const DayaBay::ReadoutHeader &readout) |
| process(const DayaBay::ReadoutPmtCrate&) process a single readout, extracting the information needed for calibration | |
| virtual StatusCode | calibrate () |
| calibrate() is called after processing many readouts. | |
| INTupleSvc * | ntupleSvc () const |
| INTupleSvc * | evtColSvc () const |
| IDataProviderSvc * | detSvc () const |
| IDataProviderSvc * | evtSvc () const |
| IIncidentSvc * | incSvc () const |
| IChronoStatSvc * | chronoSvc () const |
| IHistogramSvc * | histoSvc () const |
| IAlgContextSvc * | contextSvc () const |
| DataObject * | put (IDataProviderSvc *svc, DataObject *object, const std::string &address, const bool useRootInTES=true) const |
| DataObject * | put (DataObject *object, const std::string &address, const bool useRootInTES=true) const |
| Gaudi::Utils::GetData< TYPE >::return_type | get (IDataProviderSvc *svc, const std::string &location, const bool useRootInTES=true) const |
| Gaudi::Utils::GetData< TYPE >::return_type | get (const std::string &location, const bool useRootInTES=true) const |
| TYPE * | getDet (IDataProviderSvc *svc, const std::string &location) const |
| TYPE * | getDet (const std::string &location) const |
| bool | exist (IDataProviderSvc *svc, const std::string &location, const bool useRootInTES=true) const |
| bool | exist (const std::string &location, const bool useRootInTES=true) const |
| bool | existDet (IDataProviderSvc *svc, const std::string &location) const |
| bool | existDet (const std::string &location) const |
| TYPE * | getOrCreate (IDataProviderSvc *svc, const std::string &location, const bool useRootInTES=true) const |
| TYPE * | getOrCreate (const std::string &location, const bool useRootInTES=true) const |
| TOOL * | tool (const std::string &type, const std::string &name, const IInterface *parent=0, bool create=true) const |
| TOOL * | tool (const std::string &type, const IInterface *parent=0, bool create=true) const |
| SERVICE * | svc (const std::string &name, const bool create=true) const |
| IUpdateManagerSvc * | updMgrSvc () const |
| IDataProviderSvc * | fastContainersSvc () const |
| StatusCode | Error (const std::string &msg, const StatusCode st=StatusCode::FAILURE, const size_t mx=10) const |
| StatusCode | Warning (const std::string &msg, const StatusCode st=StatusCode::FAILURE, const size_t mx=10) const |
| StatusCode | Print (const std::string &msg, const StatusCode st=StatusCode::SUCCESS, const MSG::Level lev=MSG::INFO) const |
| StatusCode | Assert (const bool ok, const std::string &message="", const StatusCode sc=StatusCode(StatusCode::FAILURE, true)) const |
| StatusCode | Assert (const bool ok, const char *message, const StatusCode sc=StatusCode(StatusCode::FAILURE, true)) const |
| StatusCode | Exception (const std::string &msg, const GaudiException &exc, const StatusCode sc=StatusCode(StatusCode::FAILURE, true)) const |
| StatusCode | Exception (const std::string &msg, const std::exception &exc, const StatusCode sc=StatusCode(StatusCode::FAILURE, true)) const |
| StatusCode | Exception (const std::string &msg="no message", const StatusCode sc=StatusCode(StatusCode::FAILURE, true)) const |
| MsgStream & | msgStream (const MSG::Level level) const |
| MsgStream & | always () const |
| MsgStream & | fatal () const |
| MsgStream & | err () const |
| MsgStream & | error () const |
| MsgStream & | warning () const |
| MsgStream & | info () const |
| MsgStream & | debug () const |
| MsgStream & | verbose () const |
| MsgStream & | msg () const |
| const Statistics & | counters () const |
| StatEntity & | counter (const std::string &tag) const |
| MSG::Level | msgLevel () const |
| bool | msgLevel (const MSG::Level level) const |
| void | resetMsgStream () const |
| bool | typePrint () const |
| bool | propsPrint () const |
| bool | statPrint () const |
| bool | errorsPrint () const |
| long | printStat (const MSG::Level level=MSG::ALWAYS) const |
| long | printErrors (const MSG::Level level=MSG::ALWAYS) const |
| long | printProps (const MSG::Level level=MSG::ALWAYS) const |
| void | registerCondition (const std::string &condition, StatusCode(CallerClass::*mf)()=NULL) |
| void | registerCondition (const std::string &condition, CondType *&condPtrDest, StatusCode(CallerClass::*mf)()=NULL) |
| void | registerCondition (char *condition, StatusCode(CallerClass::*mf)()=NULL) |
| void | registerCondition (TargetClass *condition, StatusCode(CallerClass::*mf)()=NULL) |
| StatusCode | runUpdate () |
| TransientFastContainer< T > * | getFastContainer (const std::string &location, typename TransientFastContainer< T >::size_type initial=0) |
| StatusCode | release (const IInterface *interface) const |
| virtual unsigned long | release () |
| const std::string & | context () const |
| const std::string & | rootInTES () const |
| double | globalTimeOffset () const |
| virtual StatusCode | queryInterface (const InterfaceID &riid, void **ppvUnknown) |
| virtual unsigned long | addRef () |
| virtual const std::string & | name () const |
| virtual const std::string & | type () const |
| virtual const IInterface * | parent () const |
| virtual StatusCode | configure () |
| virtual StatusCode | start () |
| virtual StatusCode | stop () |
| virtual StatusCode | terminate () |
| virtual StatusCode | reinitialize () |
| virtual StatusCode | restart () |
| virtual Gaudi::StateMachine::State | FSMState () const |
| virtual Gaudi::StateMachine::State | targetFSMState () const |
| virtual StatusCode | sysInitialize () |
| virtual StatusCode | sysStart () |
| virtual StatusCode | sysStop () |
| virtual StatusCode | sysFinalize () |
| virtual StatusCode | sysReinitialize () |
| virtual StatusCode | sysRestart () |
| virtual StatusCode | setProperty (const Property &p) |
| virtual StatusCode | setProperty (const std::string &s) |
| virtual StatusCode | setProperty (const std::string &n, const std::string &v) |
| StatusCode | setProperty (const std::string &name, const TYPE &value) |
| virtual StatusCode | getProperty (Property *p) const |
| virtual const Property & | getProperty (const std::string &name) const |
| virtual StatusCode | getProperty (const std::string &n, std::string &v) const |
| virtual const std::vector< Property * > & | getProperties () const |
| PropertyMgr * | getPropertyMgr () |
| ISvcLocator * | serviceLocator () const |
| ISvcLocator * | svcLoc () const |
| IMessageSvc * | msgSvc () const |
| IToolSvc * | toolSvc () const |
| StatusCode | setProperties () |
| StatusCode | service (const std::string &name, T *&svc, bool createIf=true) const |
| StatusCode | service (const std::string &type, const std::string &name, T *&svc) const |
| void | declInterface (const InterfaceID &, void *) |
| Property * | declareProperty (const std::string &name, T &property, const std::string &doc="none") const |
| Property * | declareRemoteProperty (const std::string &name, IProperty *rsvc, const std::string &rname="") const |
| IAuditorSvc * | auditorSvc () const |
| IMonitorSvc * | monitorSvc () const |
| void | declareInfo (const std::string &name, const T &var, const std::string &desc) const |
| void | declareInfo (const std::string &name, const std::string &format, const void *var, int size, const std::string &desc) const |
| virtual const std::string & | type () const =0 |
| virtual const IInterface * | parent () const =0 |
| virtual StatusCode | configure ()=0 |
| virtual StatusCode | start ()=0 |
| virtual StatusCode | stop ()=0 |
| virtual StatusCode | terminate ()=0 |
| virtual StatusCode | reinitialize ()=0 |
| virtual StatusCode | restart ()=0 |
| virtual Gaudi::StateMachine::State | FSMState () const =0 |
| virtual StatusCode | sysInitialize ()=0 |
| virtual StatusCode | sysStart ()=0 |
| virtual StatusCode | sysStop ()=0 |
| virtual StatusCode | sysFinalize ()=0 |
| virtual StatusCode | sysReinitialize ()=0 |
| virtual StatusCode | sysRestart ()=0 |
| virtual unsigned long | refCount () const =0 |
| virtual const std::string & | name () const =0 |
| virtual StatusCode | queryInterface (const InterfaceID &riid, void **ppvInterface)=0 |
| virtual unsigned long | addRef ()=0 |
| virtual unsigned long | release ()=0 |
Static Public Member Functions | |
| static const InterfaceID & | interfaceID () |
| static const InterfaceID & | interfaceID () |
| static const InterfaceID & | interfaceID () |
| static const InterfaceID & | interfaceID () |
| Retrieve interface ID. | |
Public Attributes | |
| SUCCESS | |
| NO_INTERFACE | |
| VERSMISMATCH | |
| LAST_ERROR | |
Protected Types | |
| typedef std::map< std::string, StatEntity > | Statistics |
| typedef std::map< std::string, unsigned int > | Counter |
| typedef std::vector< IAlgTool * > | AlgTools |
| typedef std::pair< IInterface *, std::string > | ServiceEntry |
| typedef std::vector< ServiceEntry > | Services |
Protected Member Functions | |
| StatusCode | releaseTool (const IAlgTool *tool) const |
| StatusCode | releaseSvc (const IInterface *svc) const |
| int | outputLevel () const |
| virtual unsigned long | refCount () const |
| IntegerProperty & | outputLevelProperty () |
| void | initOutputLevel (Property &prop) |
Static Protected Attributes | |
| static const bool | IgnoreRootInTES |
| static const bool | UseRootInTES |
Private Member Functions | |
| bool | hasStats (const DayaBay::Detector &detector) |
| Check if the statistics for this detector have been initialized. | |
| StatusCode | prepareStats (const Context &context) |
| Get the current statistics for a detector, or create if needed. | |
| std::string | getPath (const DayaBay::Detector &detector) |
| std::string | getPath (const DayaBay::Detector &detector, int ring) |
| std::string | getPath (const DayaBay::FeeChannelId &channelId) |
| virtual StatusCode | calibrateRaw () |
| used for firstly calibration of PMT before calibration with cuts calibrateRaw() is called after processing many readouts. | |
Private Attributes | |
| std::string | m_cableSvcName |
| std::string | m_calibSvcName |
| std::string | m_filepath |
| ICableSvc * | m_cableSvc |
| ICalibDataSvc * | m_calibSvc |
| IStatisticsSvc * | m_statsSvc |
| std::vector< DayaBay::Detector > | m_processedDetectors |
| int | adc_buffer [max_PMT_number][max_event_number][max_FEE_hit] |
| int | tdc_buffer [max_PMT_number][max_event_number][max_FEE_hit] |
| int | adc_clock_buffer [max_PMT_number][max_event_number][max_FEE_hit] |
| double | adc_baseline_buffer [max_PMT_number][max_event_number] |
| double | m_ADC_Cut_Min |
| double | m_ADC_Cut_Max |
| std::map< int, double > | m_ADC_Cut_Min_Channel |
| std::map< int, double > | m_ADC_Cut_Max_Channel |
| double | m_TDC_Cut_Min |
| double | m_TDC_Cut_Max |
| std::map< int, double > | m_TDC_Cut_Min_Channel |
| std::map< int, double > | m_TDC_Cut_Max_Channel |
Definition at line 80 of file PmtCalibLeadingEdgeWithCuts.h.
| PmtCalibLeadingEdgeWithCuts::PmtCalibLeadingEdgeWithCuts | ( | const std::string & | type, | |
| const std::string & | name, | |||
| const IInterface * | parent | |||
| ) |
Definition at line 13 of file PmtCalibLeadingEdgeWithCuts.cc.
00016 : GaudiTool(type,name,parent) 00017 , m_cableSvc(0) 00018 , m_statsSvc(0) 00019 { 00020 declareInterface< IPmtCalibParamTool >(this) ; 00021 declareProperty("CableSvcName",m_cableSvcName="StaticCableSvc", 00022 "Name of service to map between detector, hardware, and electronic IDs"); 00023 declareProperty("CalibSvcName",m_calibSvcName="StaticCalibDataSvc", 00024 "Name of service to access FEE channel baselines"); 00025 declareProperty("FilePath",m_filepath="/file0/PmtCalibLeadingEdgeWithCuts", 00026 "File path of registered histograms."); 00027 00028 m_ADC_Cut_Min=0.25; // in ratio to single PE 00029 m_ADC_Cut_Max=-1; // in ratio to single PE 00030 00031 m_TDC_Cut_Min=15; //in ns 00032 m_TDC_Cut_Max=15; //in ns 00033 00034 for(int i=0;i<max_PMT_number;i++) 00035 for(int j=0;j<max_event_number;j++) 00036 { 00037 adc_baseline_buffer[i][j]=0; 00038 for(int k=0;k<max_FEE_hit;k++) 00039 { 00040 adc_buffer[i][j][k]=0; 00041 tdc_buffer[i][j][k]=0; 00042 adc_clock_buffer[i][j][k]=0; 00043 } 00044 } 00045 00046 }
| PmtCalibLeadingEdgeWithCuts::~PmtCalibLeadingEdgeWithCuts | ( | ) | [virtual] |
| StatusCode PmtCalibLeadingEdgeWithCuts::initialize | ( | ) | [virtual] |
Reimplemented from GaudiTool.
Definition at line 50 of file PmtCalibLeadingEdgeWithCuts.cc.
00051 { 00052 info() << "initialize()" << endreq; 00053 StatusCode sc = this->GaudiTool::initialize(); 00054 if( sc != StatusCode::SUCCESS ) return sc; 00055 00056 // Get Cable Service 00057 m_cableSvc = svc<ICableSvc>(m_cableSvcName,true); 00058 if(!m_cableSvc){ 00059 error() << "Failed to access cable svc: " << m_cableSvcName << endreq; 00060 return StatusCode::FAILURE; 00061 } 00062 00063 // Get Calibration Service 00064 m_calibSvc = svc<ICalibDataSvc>(m_calibSvcName,true); 00065 if(!m_calibSvc){ 00066 error() << "Failed to access calibration svc: " << m_calibSvcName << endreq; 00067 return StatusCode::FAILURE; 00068 } 00069 00070 // Get the histogram service 00071 m_statsSvc = svc<IStatisticsSvc>("StatisticsSvc",true); 00072 if(!m_statsSvc) { 00073 error() << " No StatisticsSvc available." << endreq; 00074 return StatusCode::FAILURE; 00075 } 00076 00077 // Initialize data from input file. 00078 // Processed detectors are identified by directory name. 00079 std::vector<std::string> detDirs = m_statsSvc->getSubFolders(m_filepath); 00080 std::vector<std::string>::iterator dirIter, dirEnd = detDirs.end(); 00081 for(dirIter=detDirs.begin(); dirIter != dirEnd; dirIter++){ 00082 info() << "Processing input directory: " << m_filepath << "/" 00083 << *dirIter << endreq; 00084 DayaBay::Detector detector( 00085 DayaBay::Detector::siteDetPackedFromString(*dirIter) ); 00086 if( detector.site() == Site::kUnknown 00087 || detector.detectorId() == DetectorId::kUnknown ){ 00088 warning() << "Input directory: " << m_filepath << "/" 00089 << *dirIter << " not processed." << endreq; 00090 } 00091 m_processedDetectors.push_back(detector); 00092 } 00093 00094 return StatusCode::SUCCESS; 00095 }
| StatusCode PmtCalibLeadingEdgeWithCuts::finalize | ( | ) | [virtual] |
Reimplemented from GaudiTool.
Definition at line 97 of file PmtCalibLeadingEdgeWithCuts.cc.
00098 { 00099 StatusCode sc = this->GaudiTool::finalize(); 00100 if( sc != StatusCode::SUCCESS ) return sc; 00101 00102 return StatusCode::SUCCESS; 00103 }
| StatusCode PmtCalibLeadingEdgeWithCuts::process | ( | const DayaBay::ReadoutHeader & | readout | ) | [virtual] |
process(const DayaBay::ReadoutPmtCrate&) process a single readout, extracting the information needed for calibration
Implements IPmtCalibParamTool.
Definition at line 107 of file PmtCalibLeadingEdgeWithCuts.cc.
00108 { 00109 /* 00110 DayaBay::Detector detector = readout.detector(); 00111 // Only process AD detectors 00112 if(detector.detectorId() != DetectorId::kAD1 00113 && detector.detectorId() != DetectorId::kAD2 00114 && detector.detectorId() != DetectorId::kAD3 00115 && detector.detectorId() != DetectorId::kAD4){ 00116 debug() << "process(): Do not process detector: " << detector.detName() 00117 << endreq; 00118 return StatusCode::SUCCESS; 00119 } 00120 00121 const DayaBay::ReadoutHeader* header = readout.header(); 00122 if(!header){ 00123 error() << "process(): Readout has no header!" << endreq; 00124 return StatusCode::FAILURE; 00125 } 00126 Context context = header->context(); 00127 if( !hasStats(detector) ){ 00128 // New detector; initialize all the detector statistics 00129 StatusCode sc = prepareStats(context); 00130 if( sc != StatusCode::SUCCESS ) return sc; 00131 } 00132 ServiceMode svcMode(context, 0); 00133 00134 // Retrieve the parameters and histograms 00135 TObject* nReadoutsObj = 00136 m_statsSvc->get(this->getPath(detector) + "/nReadouts"); 00137 TParameter<int>* nReadoutsPar= dynamic_cast<TParameter<int>*>(nReadoutsObj); 00138 if(!nReadoutsPar) return StatusCode::FAILURE; 00139 TH1F* adcSumH = m_statsSvc->getTH1F( this->getPath(detector) + "/adcSum"); 00140 if(!adcSumH) return StatusCode::FAILURE; 00141 TH1F* tdcMedianH = m_statsSvc->getTH1F( this->getPath(detector)+"/tdcMedian"); 00142 if(!tdcMedianH) return StatusCode::FAILURE; 00143 00144 // Add the current readout to the calibration statistics 00145 std::vector<int> readoutTdc; 00146 double adcSum = 0; 00147 int nReadouts = nReadoutsPar->GetVal(); 00148 // Loop over channels 00149 DayaBay::ReadoutPmtCrate::PmtChannelReadouts::const_iterator chanIter, 00150 chanEnd = readout.channelReadout().end(); 00151 for(chanIter = readout.channelReadout().begin(); chanIter != chanEnd; 00152 chanIter++){ 00153 const DayaBay::ReadoutPmtChannel& pmtChan = chanIter->second; 00154 const DayaBay::FeeChannelId& chanId = pmtChan.channelId(); 00155 const DayaBay::FeeCalibData* feeCalib = m_calibSvc->feeCalibData(chanId, 00156 svcMode); 00157 00158 // Analyze the channel data 00159 // Determine ring and column 00160 DayaBay::AdPmtSensor pmtId = m_cableSvc->adPmtSensor(chanId, svcMode); 00161 int ring = pmtId.ring(); 00162 int column = pmtId.column(); 00163 int chan_num=(ring*24+column)%max_PMT_number; 00164 00165 //info()<<"test 7 chan_num="<<chan_num<<", ring="<<ring<<",column="<<column<<endreq; 00166 00167 00168 // Get Parameters and Histograms 00169 TObject* nHitsObj = m_statsSvc->get(this->getPath(chanId) + "/nHits"); 00170 TParameter<int>* nHitsPar = dynamic_cast<TParameter<int>*>(nHitsObj); 00171 if(!nHitsPar) return StatusCode::FAILURE; 00172 TH1F* tdcRawH = m_statsSvc->getTH1F( this->getPath(chanId) + "/tdcRaw"); 00173 if(!tdcRawH) return StatusCode::FAILURE; 00174 TH1F* adcRawH = m_statsSvc->getTH1F( this->getPath(chanId) + "/adcRaw"); 00175 if(!adcRawH) return StatusCode::FAILURE; 00176 TH1F* adcH = m_statsSvc->getTH1F( this->getPath(chanId) + "/adc"); 00177 if(!adcH) return StatusCode::FAILURE; 00178 TH1F* adcByClockH = m_statsSvc->getTH1F( this->getPath(chanId) 00179 + "/adcByClock"); 00180 if(!adcByClockH) return StatusCode::FAILURE; 00181 00182 std::vector<int>::const_iterator tdcIter, tdcEnd = pmtChan.tdc().end(); 00183 // Loop over tdc values 00184 int tdc_hit_count=0; 00185 for(tdcIter = pmtChan.tdc().begin(); tdcIter!=tdcEnd; tdcIter++){ 00186 int tdc = *tdcIter; 00187 readoutTdc.push_back(tdc); 00188 tdcRawH->Fill(tdc); 00189 tdc_buffer[chan_num][nReadouts%max_event_number][tdc_hit_count%max_FEE_hit]=tdc; 00190 tdc_hit_count++; 00191 } 00192 // Loop over adc values 00193 double adcBaseline = feeCalib->m_adcBaselineLow; 00194 adc_baseline_buffer[chan_num][nReadouts%max_event_number]=adcBaseline; 00195 int adc_hit_count=0; 00196 for(int i=0; i<pmtChan.size(); i++) { 00197 int adcClock = pmtChan.adcCycle(i); 00198 int adc = pmtChan.adc(i); 00199 adcRawH->Fill(adc); 00200 adcH->Fill(adc-adcBaseline); 00201 adcByClockH->Fill(adcClock,adc-adcBaseline); 00202 adcSum += adc-adcBaseline; 00203 adc_buffer[chan_num][nReadouts%max_event_number][adc_hit_count%max_FEE_hit]=adc; 00204 //info()<<"test 4 chan_num="<<chan_num<<",nReadouts= "<<nReadouts<<", adc_hit_count="<<adc_hit_count<<",adc="<<adc_buffer[chan_num][nReadouts%max_event_number][adc_hit_count%max_FEE_hit]<<endreq; 00205 adc_clock_buffer[chan_num][nReadouts%max_event_number][adc_hit_count%max_FEE_hit]=adcClock; 00206 adc_hit_count++; 00207 } 00208 // Increment number of hits on this channel 00209 nHitsPar->SetVal( nHitsPar->GetVal() + 1 ); 00210 } 00211 00212 00213 if(readoutTdc.size() > 0){ 00214 // Find the median TDC 00215 std::sort(readoutTdc.begin(), readoutTdc.end()); 00216 int medianIdx = readoutTdc.size()/2; 00217 int medianTdc = readoutTdc[medianIdx]; 00218 for(chanIter = readout.channelReadout().begin(); chanIter != chanEnd; 00219 chanIter++){ 00220 const DayaBay::ReadoutPmtChannel& pmtChan = chanIter->second; 00221 const DayaBay::FeeChannelId& chanId = pmtChan.channelId(); 00222 TH1F* tdcByMedianH = 00223 m_statsSvc->getTH1F( this->getPath(chanId) + "/tdcByMedian"); 00224 if(!tdcByMedianH) return StatusCode::FAILURE; 00225 std::vector<int>::const_iterator tdcIter, tdcEnd = pmtChan.tdc().end(); 00226 // Loop over tdc values 00227 for(tdcIter = pmtChan.tdc().begin(); tdcIter!=tdcEnd; tdcIter++){ 00228 int tdc = *tdcIter; 00229 tdcByMedianH->Fill(tdc-medianTdc); 00230 } 00231 } 00232 tdcMedianH->Fill(medianTdc); 00233 } 00234 00235 adcSumH->Fill(adcSum); 00236 // Increment number of readouts from this detector 00237 nReadoutsPar->SetVal( nReadoutsPar->GetVal() + 1 ); 00238 00239 00240 00241 00242 //info() << "test first data read over" << endreq; 00243 00244 // used for preliminary calibrate the PMT and get preliminary results 00245 // the condition for this is the nReadouts+1==max_event_number 00246 // if the total nReadouts+1 of whole analysis is less than max_event_number, the calibrateRaw() will invoke with finalize() 00247 // wangzhm@ihep.ac.cn 2009/11/29 00248 if(nReadouts+1==max_event_number) 00249 if(calibrateRaw()==StatusCode::FAILURE) 00250 { 00251 warning() << "calibrateRaw() error" 00252 << endreq; 00253 return StatusCode::FAILURE; 00254 } 00255 00256 //info() << "test first calibration over" << endreq; 00257 00258 00260 // read the data again with the cuts from firstly calibration 00261 // wangzhm@ihep.ac.cn 2009/11/28 00263 00264 00265 // the condition for this is the nReadouts>max_event_number 00266 // if the total nReadouts+1 of whole analysis is less than max_event_number, the calibrateRaw() will invoke with finalize() 00267 // wangzhm@ihep.ac.cn 2009/11/29 00268 if(nReadouts+1>max_event_number) 00269 { 00270 info()<<"test 1 nReadouts="<<nReadouts<<endreq; 00271 00272 TH1F* adcSumHWithCuts = m_statsSvc->getTH1F( this->getPath(detector) + "/adcSumWithCuts"); 00273 if(!adcSumHWithCuts) return StatusCode::FAILURE; 00274 TH1F* tdcMedianHWithCuts = m_statsSvc->getTH1F( this->getPath(detector)+"/tdcMedianWithCuts"); 00275 if(!tdcMedianHWithCuts) return StatusCode::FAILURE; 00276 00277 // Add the current readout to the calibration statistics 00278 std::vector<int> readoutTdcWithCuts; 00279 double adcSumWithCuts = 0; 00280 // Loop over channels 00281 for(chanIter = readout.channelReadout().begin(); chanIter != chanEnd; 00282 chanIter++){ 00283 const DayaBay::ReadoutPmtChannel& pmtChan = chanIter->second; 00284 const DayaBay::FeeChannelId& chanId = pmtChan.channelId(); 00285 const DayaBay::FeeCalibData* feeCalib = m_calibSvc->feeCalibData(chanId, 00286 svcMode); 00287 00288 // Analyze the channel data 00289 // Determine ring and column 00290 DayaBay::AdPmtSensor pmtId = m_cableSvc->adPmtSensor(chanId, svcMode); 00291 int ring = pmtId.ring(); 00292 int column = pmtId.column(); 00293 int channel_number=ring*100+column; 00294 00295 //info()<<"process ring="<<ring<<" column="<<column<<endreq; 00296 00297 00298 // Get Parameters and Histograms 00299 TObject* nHitsObjWithCuts = m_statsSvc->get(this->getPath(chanId) + "/nHitsWithCuts"); 00300 TParameter<int>* nHitsParWithCuts = dynamic_cast<TParameter<int>*>(nHitsObjWithCuts); 00301 if(!nHitsParWithCuts) return StatusCode::FAILURE; 00302 TH1F* tdcHWithCuts = m_statsSvc->getTH1F( this->getPath(chanId) + "/tdcWithCuts"); 00303 if(!tdcHWithCuts) return StatusCode::FAILURE; 00304 TH1F* adcRawHWithCuts = m_statsSvc->getTH1F( this->getPath(chanId) + "/adcRawWithCuts"); 00305 if(!adcRawHWithCuts) return StatusCode::FAILURE; 00306 TH1F* adcHWithCuts = m_statsSvc->getTH1F( this->getPath(chanId) + "/adcWithCuts"); 00307 if(!adcHWithCuts) return StatusCode::FAILURE; 00308 TH1F* adcByClockHWithCuts = m_statsSvc->getTH1F( this->getPath(chanId) 00309 + "/adcByClockWithCuts"); 00310 if(!adcByClockHWithCuts) return StatusCode::FAILURE; 00311 00312 bool hit_ok_flag=0; 00313 double adcBaseline = feeCalib->m_adcBaselineLow; 00314 // Loop over tdc values 00315 for(int i=0; i<pmtChan.size(); i++) { 00316 int tdc = pmtChan.tdc(i); 00317 00318 int adcClock = pmtChan.adcCycle(i); 00319 int adc = pmtChan.adc(i); 00320 00321 if(tdc>=m_TDC_Cut_Min_Channel[channel_number]&&tdc<=m_TDC_Cut_Max_Channel[channel_number]) 00322 if(adc>=m_ADC_Cut_Min_Channel[channel_number]&&adc<=m_ADC_Cut_Max_Channel[channel_number]) 00323 { 00324 readoutTdcWithCuts.push_back(tdc); 00325 tdcHWithCuts->Fill(tdc); 00326 00327 adcRawHWithCuts->Fill(adc); 00328 adcHWithCuts->Fill(adc-adcBaseline); 00329 adcByClockHWithCuts->Fill(adcClock,adc-adcBaseline); 00330 adcSumWithCuts += adc-adcBaseline; 00331 hit_ok_flag=1; 00332 } 00333 00334 } 00335 // Increment number of hits on this channel 00336 if(hit_ok_flag) 00337 nHitsParWithCuts->SetVal( nHitsParWithCuts->GetVal() + 1 ); 00338 } 00339 00340 if(readoutTdcWithCuts.size() > 0){ 00341 // Find the median TDC 00342 std::sort(readoutTdcWithCuts.begin(), readoutTdcWithCuts.end()); 00343 int medianIdx = readoutTdcWithCuts.size()/2; 00344 int medianTdc = readoutTdcWithCuts[medianIdx]; 00345 for(chanIter = readout.channelReadout().begin(); chanIter != chanEnd; 00346 chanIter++){ 00347 const DayaBay::ReadoutPmtChannel& pmtChan = chanIter->second; 00348 const DayaBay::FeeChannelId& chanId = pmtChan.channelId(); 00349 00350 // Analyze the channel data 00351 // Determine ring and column 00352 DayaBay::AdPmtSensor pmtId = m_cableSvc->adPmtSensor(chanId, svcMode); 00353 int ring = pmtId.ring(); 00354 int column = pmtId.column(); 00355 int channel_number=ring*100+column; 00356 00357 00358 TH1F* tdcByMedianHWithCuts = 00359 m_statsSvc->getTH1F( this->getPath(chanId) + "/tdcByMedianWithCuts"); 00360 if(!tdcByMedianHWithCuts) return StatusCode::FAILURE; 00361 // Loop over tdc values 00362 for(int i=0; i<pmtChan.size(); i++) { 00363 int tdc = pmtChan.tdc(i); 00364 int adc = pmtChan.adc(i); 00365 00366 if(tdc>=m_TDC_Cut_Min_Channel[channel_number]&&tdc<=m_TDC_Cut_Max_Channel[channel_number]) 00367 if(adc>=m_ADC_Cut_Min_Channel[channel_number]&&adc<=m_ADC_Cut_Max_Channel[channel_number]) 00368 tdcByMedianHWithCuts->Fill(tdc-medianTdc); 00369 } 00370 } 00371 tdcMedianHWithCuts->Fill(medianTdc); 00372 } 00373 00374 adcSumHWithCuts->Fill(adcSumWithCuts); 00375 00376 00377 }//over for if(nReadouts>=max_event_number) 00378 //info() << "test second read data over" << endreq; 00379 00382 00383 */ 00384 return StatusCode::SUCCESS; 00385 }
| StatusCode PmtCalibLeadingEdgeWithCuts::calibrate | ( | ) | [virtual] |
calibrate() is called after processing many readouts.
This method is responsible for calculating the calibration parameters.
Implements IPmtCalibParamTool.
Definition at line 845 of file PmtCalibLeadingEdgeWithCuts.cc.
00845 { 00846 // Loop over detectors in summary data, and calculate parameters 00847 00848 // before the final calibrate of PMT 00849 // we need to fit the raw data without any cuts 00850 // and if the total nreadouts number less than max_event_number, 00851 // the data with cuts need to be fill firstly which will realized in the same program 00852 // 2009/11/29 00853 if(calibrateRaw()==StatusCode::FAILURE) 00854 { 00855 warning() << "calibrateRaw() error" 00856 << endreq; 00857 return StatusCode::FAILURE; 00858 } 00859 00860 00861 std::vector<DayaBay::Detector>::iterator detIter, 00862 detEnd = m_processedDetectors.end(); 00863 for(detIter = m_processedDetectors.begin(); detIter != detEnd; detIter++){ 00864 DayaBay::Detector detector = *detIter; 00865 00866 // Retrieve the data description and histograms 00867 TObject* nReadoutsObj = 00868 m_statsSvc->get(this->getPath(detector) + "/nReadouts"); 00869 TParameter<int>* nReadoutsPar= dynamic_cast<TParameter<int>*>(nReadoutsObj); 00870 if(!nReadoutsPar) return StatusCode::FAILURE; 00871 00872 TH1F* meanOccupancyHWithCuts = m_statsSvc->getTH1F( this->getPath(detector) 00873 +"/meanOccupancyWithCuts"); 00874 if(!meanOccupancyHWithCuts) return StatusCode::FAILURE; 00875 00876 TH1F* meanAdcHWithCuts = m_statsSvc->getTH1F( this->getPath(detector) 00877 +"/meanAdcWithCuts"); 00878 if(!meanAdcHWithCuts) return StatusCode::FAILURE; 00879 00880 00881 TH1F* gainChannelHWithCuts = m_statsSvc->getTH1F( this->getPath(detector) 00882 +"/gainChannelWithCuts"); 00883 if(!gainChannelHWithCuts) return StatusCode::FAILURE; 00884 00885 TH1F* tdcOffsetChannelHWithCuts = m_statsSvc->getTH1F( this->getPath(detector) 00886 +"/tdcOffsetChannelWithCuts"); 00887 if(!tdcOffsetChannelHWithCuts) return StatusCode::FAILURE; 00888 TH1F* occupancyChannelHWithCuts = m_statsSvc->getTH1F( this->getPath(detector) 00889 +"/occupancyChannelWithCuts"); 00890 if(!occupancyChannelHWithCuts) return StatusCode::FAILURE; 00891 00892 00893 TH1F* meanTdcOffsetHWithCuts = m_statsSvc->getTH1F( this->getPath(detector) 00894 +"/meanTdcOffsetWithCuts"); 00895 if(!meanTdcOffsetHWithCuts) return StatusCode::FAILURE; 00896 00897 TObject* simFlagObj = m_statsSvc->get(this->getPath(detector) +"/simFlag"); 00898 TParameter<int>* simFlagPar = dynamic_cast<TParameter<int>*>(simFlagObj); 00899 if(!simFlagPar) return StatusCode::FAILURE; 00900 00901 TObject* timeSecObj = m_statsSvc->get(this->getPath(detector) +"/timeSec"); 00902 TParameter<int>* timeSecPar = dynamic_cast<TParameter<int>*>(timeSecObj); 00903 if(!timeSecPar) return StatusCode::FAILURE; 00904 00905 TObject* timeNanoSecObj = 00906 m_statsSvc->get(this->getPath(detector) +"/timeNanoSec"); 00907 TParameter<int>* timeNanoSecPar = 00908 dynamic_cast<TParameter<int>*>(timeNanoSecObj); 00909 if(!timeNanoSecPar) return StatusCode::FAILURE; 00910 00911 Site::Site_t site = detector.site(); 00912 DetectorId::DetectorId_t detId = detector.detectorId(); 00913 SimFlag::SimFlag_t simFlag = (SimFlag::SimFlag_t)(simFlagPar->GetVal()); 00914 time_t timeSec = (time_t)(timeSecPar->GetVal()); 00915 int timeNanoSec = timeNanoSecPar->GetVal(); 00916 int nReadouts = nReadoutsPar->GetVal(); 00917 00918 // Context, ServiceMode 00919 Context context(site,simFlag,TimeStamp(timeSec,timeNanoSec),detId); 00920 ServiceMode svcMode(context, 0); 00921 00922 // Prepare data for ring-to-ring comparison 00923 int nRings = 8; 00924 std::map<int,double> meanOccRing; 00925 std::map<int,double> meanOccWeight; 00926 std::map<int,double> meanTdcRing; 00927 std::map<int,double> meanTdcWeight; 00928 for(int ring = 1; ring <= nRings; ring++){ 00929 meanOccRing[ring] = 0.0; 00930 meanOccWeight[ring] = 0.0; 00931 meanTdcRing[ring] = 0.0; 00932 meanTdcWeight[ring] = 0.0; 00933 } 00934 double minValidOcc = 0.01; // Cut on PMTs with bad/low occupancy 00935 std::vector<DayaBay::FeeChannelId> badChannels; 00936 00937 // Get list of all FEE channels 00938 const std::vector<DayaBay::FeeChannelId>& channelList 00939 = m_cableSvc->feeChannelIds( svcMode ); 00940 std::vector<DayaBay::FeeChannelId>::const_iterator chanIter, 00941 chanEnd = channelList.end(); 00942 // Initialize statistics for each channel 00943 for(chanIter = channelList.begin(); chanIter != chanEnd; chanIter++){ 00944 DayaBay::FeeChannelId chanId = *chanIter; 00945 00946 // Analyze the channel data 00947 // Determine ring and column 00948 DayaBay::AdPmtSensor pmtId = 00949 m_cableSvc->adPmtSensor(chanId, svcMode); 00950 int ring = pmtId.ring(); 00951 int column = pmtId.column(); 00952 00953 // Channel status 00954 bool occupancyOk = true; 00955 bool gainOk = true; 00956 00957 // Get Parameters and Histograms 00958 TObject* nHitsObjWithCuts = m_statsSvc->get(this->getPath(chanId) + "/nHitsWithCuts"); 00959 TParameter<int>* nHitsParWithCuts = dynamic_cast<TParameter<int>*>(nHitsObjWithCuts); 00960 if(!nHitsParWithCuts) return StatusCode::FAILURE; 00961 TH1F* tdcByMedianHWithCuts= m_statsSvc->getTH1F( this->getPath(chanId) 00962 +"/tdcByMedianWithCuts"); 00963 if(!tdcByMedianHWithCuts) return StatusCode::FAILURE; 00964 TH1F* adcHWithCuts = m_statsSvc->getTH1F( this->getPath(chanId) + "/adcWithCuts"); 00965 if(!adcHWithCuts) return StatusCode::FAILURE; 00966 TH1F* occupancyHWithCuts = m_statsSvc->getTH1F( this->getPath(detector,ring) 00967 +"/occupancyWithCuts"); 00968 if(!occupancyHWithCuts) return StatusCode::FAILURE; 00969 TH1F* tdcOffsetHWithCuts = m_statsSvc->getTH1F( this->getPath(detector,ring) 00970 +"/tdcOffsetWithCuts"); 00971 if(!tdcOffsetHWithCuts) return StatusCode::FAILURE; 00972 TH1F* adcMedianHWithCuts = m_statsSvc->getTH1F( this->getPath(detector,ring) 00973 +"/adcMedianWithCuts"); 00974 if(!adcMedianHWithCuts) return StatusCode::FAILURE; 00975 TH1F* adcSigmaHWithCuts = m_statsSvc->getTH1F( this->getPath(detector,ring) 00976 +"/adcSigmaWithCuts"); 00977 if(!adcSigmaHWithCuts) return StatusCode::FAILURE; 00978 00979 // Channel Occupancy 00980 double occupancy = 0; 00981 double occupancyUncert = 0; 00982 { 00983 int nHits = nHitsParWithCuts->GetVal(); 00984 if(nReadouts > 0){ 00985 occupancy = nHits / ((double) nReadouts); 00986 occupancyUncert = 00987 sqrt( occupancy*(1-occupancy)/((double) nReadouts) ); 00988 } 00989 occupancyHWithCuts->SetBinContent(occupancyHWithCuts->FindBin(column),occupancy); 00990 occupancyHWithCuts->SetBinError(occupancyHWithCuts->FindBin(column),occupancyUncert); 00991 if(occupancy < minValidOcc) occupancyOk = false; 00992 00993 occupancyChannelHWithCuts->SetBinContent(ring*24+column,occupancy); 00994 occupancyChannelHWithCuts->SetBinError(ring*24+column,occupancyUncert); 00995 00996 00997 } 00998 // TDC offset 00999 double tdcOffset = 0; 01000 double tdcOffsetUncert = 0; 01001 { 01002 // Use a fit to the leading edge to determine the time offset 01003 int maxBin = tdcByMedianHWithCuts->GetMaximumBin(); 01004 double fitMin = tdcByMedianHWithCuts->GetBinCenter(maxBin - 2); 01005 double fitMax = tdcByMedianHWithCuts->GetBinCenter(maxBin + 10); 01006 tdcByMedianHWithCuts->Fit("gaus","","",fitMin, fitMax); 01007 TF1* fitFunc = tdcByMedianHWithCuts->GetFunction("gaus"); 01008 if( fitFunc ){ 01009 tdcOffset = fitFunc->GetParameter(1); 01010 tdcOffsetUncert = fitFunc->GetParError(1); 01011 01012 if( fitFunc->GetNDF() < 3 ) 01013 // Catch bad fits 01014 tdcOffsetUncert = 0; 01015 else 01016 // Scale uncertainty by the quality of the fit 01017 tdcOffsetUncert *= sqrt(fitFunc->GetChisquare() / fitFunc->GetNDF()); 01018 01019 tdcOffsetHWithCuts->SetBinContent(tdcOffsetHWithCuts->FindBin(column),tdcOffset); 01020 tdcOffsetHWithCuts->SetBinError(tdcOffsetHWithCuts->FindBin(column),tdcOffsetUncert); 01021 01022 tdcOffsetChannelHWithCuts->SetBinContent(ring*24+column,tdcOffset); 01023 tdcOffsetChannelHWithCuts->SetBinError(ring*24+column,tdcOffsetUncert); 01024 01025 01026 }else{ 01027 warning() << "PMT " << pmtId << " has no TDC data for fitting." 01028 << endreq; 01029 } 01030 } 01031 01032 // ADC median and sigma 01033 double adcMean = 0; 01034 double adcMeanUncert = 0; 01035 double adcSigma = 0; 01036 double adcSigmaUncert = 0; 01037 { 01038 // Find ADC SPE peak, width with simple gaussian 01039 TF1 adcFit("adcFit","gaus"); 01040 double adcPeak = adcHWithCuts->GetBinCenter(adcHWithCuts->GetMaximumBin()); 01041 double adcSig = 5.0; 01042 adcHWithCuts->Fit(&adcFit,"","",adcPeak - adcSig, adcPeak + adcSig); 01043 adcMean = adcFit.GetParameter(1); 01044 if(fabs(adcMean-adcPeak) > adcSig){ 01045 warning() << "Issues with fit on ring, column: " << ring << " / " 01046 << column << " mean/peak: " << adcMean << " / " << adcPeak 01047 << endreq; 01048 gainOk = false; 01049 } 01050 adcMeanUncert = adcFit.GetParError(1); 01051 adcSigma = adcFit.GetParameter(2); 01052 adcSigmaUncert = adcFit.GetParError(2); 01053 adcMedianHWithCuts->SetBinContent(adcMedianHWithCuts->FindBin(column),adcMean); 01054 adcMedianHWithCuts->SetBinError(adcMedianHWithCuts->FindBin(column),adcMeanUncert); 01055 adcSigmaHWithCuts->SetBinContent(adcSigmaHWithCuts->FindBin(column),adcSigma); 01056 adcSigmaHWithCuts->SetBinError(adcSigmaHWithCuts->FindBin(column),adcSigmaUncert); 01057 01058 meanAdcHWithCuts->SetBinContent(ring,adcMean); 01059 meanAdcHWithCuts->SetBinError(ring, adcMeanUncert); 01060 01061 gainChannelHWithCuts->SetBinContent(ring*24+column,adcMean); 01062 gainChannelHWithCuts->SetBinError(ring*24+column,adcMeanUncert); 01063 01064 } 01065 01066 // Record ring averages, ignoring bad channels 01067 if( occupancyOk && gainOk ){ 01068 if( occupancyUncert > 0 ){ 01069 meanOccRing[ring] += occupancy / (occupancyUncert*occupancyUncert); 01070 meanOccWeight[ring] += 1.0 / (occupancyUncert*occupancyUncert); 01071 } 01072 if( tdcOffsetUncert > 0 ){ 01073 meanTdcRing[ring] += tdcOffset / (tdcOffsetUncert*tdcOffsetUncert); 01074 meanTdcWeight[ring] += 1.0 / (tdcOffsetUncert*tdcOffsetUncert); 01075 } 01076 } 01077 01078 } 01079 01080 // Record mean time offset and mean occupancy by ring 01081 for(int ring = 1; ring <= nRings; ring++){ 01082 if( meanOccWeight[ring] > 0 ){ 01083 meanOccupancyHWithCuts->SetBinContent(ring, 01084 meanOccRing[ring] / meanOccWeight[ring]); 01085 meanOccupancyHWithCuts->SetBinError(ring, 01086 sqrt( 1.0 / meanOccWeight[ring] )); 01087 } 01088 if( meanTdcWeight[ring] > 0 ){ 01089 meanTdcOffsetHWithCuts->SetBinContent(ring, 01090 meanTdcRing[ring] / meanTdcWeight[ring]); 01091 meanTdcOffsetHWithCuts->SetBinError(ring, 01092 sqrt( 1.0 / meanTdcWeight[ring] )); 01093 } 01094 } // Loop over rings 01095 01096 } // Loop over detectors 01097 01098 info() << "test second calibrate over" << endreq; 01099 01100 return StatusCode::SUCCESS; 01101 01102 }
| bool PmtCalibLeadingEdgeWithCuts::hasStats | ( | const DayaBay::Detector & | detector | ) | [private] |
Check if the statistics for this detector have been initialized.
Definition at line 1108 of file PmtCalibLeadingEdgeWithCuts.cc.
01108 { 01109 // Check if statistics have been initialized for this detector 01110 return (std::find(m_processedDetectors.begin(), 01111 m_processedDetectors.end(), 01112 detector) 01113 != m_processedDetectors.end()); 01114 }
| StatusCode PmtCalibLeadingEdgeWithCuts::prepareStats | ( | const Context & | context | ) | [private] |
Get the current statistics for a detector, or create if needed.
Definition at line 1116 of file PmtCalibLeadingEdgeWithCuts.cc.
01116 { 01117 // Create the histograms and parameters for this detector's statistics 01118 DayaBay::Detector detector(context.GetSite(), context.GetDetId()); 01119 01120 // Site 01121 { 01122 std::string name = "site"; 01123 std::ostringstream path; 01124 path << this->getPath(detector) << "/" << name; 01125 TParameter<int>* par = new TParameter<int>(name.c_str(), 01126 context.GetSite()); 01127 if( m_statsSvc->put(path.str(),par).isFailure() ) { 01128 error() << "prepareStats(): Could not register " << name 01129 << " at " << path << endreq; 01130 delete par; 01131 par = 0; 01132 return StatusCode::FAILURE; 01133 } 01134 } 01135 01136 // Detector ID 01137 { 01138 std::string name = "detectorId"; 01139 std::ostringstream path; 01140 path << this->getPath(detector) << "/" << name; 01141 TParameter<int>* par = new TParameter<int>(name.c_str(), 01142 context.GetDetId()); 01143 if( m_statsSvc->put(path.str(),par).isFailure() ) { 01144 error() << "prepareStats(): Could not register " << name 01145 << " at " << path << endreq; 01146 delete par; 01147 par = 0; 01148 return StatusCode::FAILURE; 01149 } 01150 } 01151 01152 // SimFlag 01153 { 01154 std::string name = "simFlag"; 01155 std::ostringstream path; 01156 path << this->getPath(detector) << "/" << name; 01157 TParameter<int>* par = new TParameter<int>(name.c_str(), 01158 context.GetSimFlag()); 01159 if( m_statsSvc->put(path.str(),par).isFailure() ) { 01160 error() << "prepareStats(): Could not register " << name 01161 << " at " << path << endreq; 01162 delete par; 01163 par = 0; 01164 return StatusCode::FAILURE; 01165 } 01166 } 01167 01168 // Timestamp: seconds 01169 { 01170 std::string name = "timeSec"; 01171 std::ostringstream path; 01172 path << this->getPath(detector) << "/" << name; 01173 TParameter<int>* par = new TParameter<int>(name.c_str(), 01174 context.GetTimeStamp().GetSec()); 01175 if( m_statsSvc->put(path.str(),par).isFailure() ) { 01176 error() << "prepareStats(): Could not register " << name 01177 << " at " << path << endreq; 01178 delete par; 01179 par = 0; 01180 return StatusCode::FAILURE; 01181 } 01182 } 01183 01184 // Timestamp: nanoseconds 01185 { 01186 std::string name = "timeNanoSec"; 01187 std::ostringstream path; 01188 path << this->getPath(detector) << "/" << name; 01189 TParameter<int>* par = 01190 new TParameter<int>(name.c_str(), context.GetTimeStamp().GetNanoSec()); 01191 if( m_statsSvc->put(path.str(),par).isFailure() ) { 01192 error() << "prepareStats(): Could not register " << name 01193 << " at " << path << endreq; 01194 delete par; 01195 par = 0; 01196 return StatusCode::FAILURE; 01197 } 01198 } 01199 01200 // Number of Readouts 01201 { 01202 std::string name = "nReadouts"; 01203 std::ostringstream path; 01204 path << this->getPath(detector) << "/" << name; 01205 TParameter<int>* par = new TParameter<int>(name.c_str(), 0); 01206 if( m_statsSvc->put(path.str(),par).isFailure() ) { 01207 error() << "prepareStats(): Could not register " << name 01208 << " at " << path << endreq; 01209 delete par; 01210 par = 0; 01211 return StatusCode::FAILURE; 01212 } 01213 } 01214 01215 // ADC Sum spectrum 01216 { 01217 std::ostringstream title, path; 01218 std::string name = "adcSum"; 01219 title << "ADC Sum for readouts in " << detector.detName(); 01220 path << this->getPath(detector) << "/" << name; 01221 TH1F* adcSum = new TH1F(name.c_str(),title.str().c_str(), 01222 2000,0,20000); 01223 adcSum->GetXaxis()->SetTitle("ADC Sum value"); 01224 adcSum->GetYaxis()->SetTitle("Entries"); 01225 // DEBUG 01226 info() << "name= " << adcSum->GetName() 01227 << " at path = \"" << path.str() << "\"" << endreq; 01228 if( m_statsSvc->put(path.str(),adcSum).isFailure() ) { 01229 error() << "prepareStats(): Could not register " << path << endreq; 01230 delete adcSum; 01231 return StatusCode::FAILURE; 01232 } 01233 } 01234 01235 // ADC Sum spectrum with cuts 01236 { 01237 std::ostringstream title, path; 01238 std::string name = "adcSumWithCuts"; 01239 title << "ADC Sum for readouts in " << detector.detName(); 01240 path << this->getPath(detector) << "/" << name; 01241 TH1F* adcSumWithCuts = new TH1F(name.c_str(),title.str().c_str(), 01242 2000,0,20000); 01243 adcSumWithCuts->GetXaxis()->SetTitle("ADC Sum value With Cuts"); 01244 adcSumWithCuts->GetYaxis()->SetTitle("Entries"); 01245 // DEBUG 01246 info() << "name= " << adcSumWithCuts->GetName() 01247 << " at path = \"" << path.str() << "\"" << endreq; 01248 if( m_statsSvc->put(path.str(),adcSumWithCuts).isFailure() ) { 01249 error() << "prepareStats(): Could not register " << path << endreq; 01250 delete adcSumWithCuts; 01251 return StatusCode::FAILURE; 01252 } 01253 } 01254 01255 01256 // Gain vs channel spectrum 01257 { 01258 std::ostringstream title, path; 01259 std::string name = "gainChannel"; 01260 title << "gainChannel for readouts in " << detector.detName(); 01261 path << this->getPath(detector) << "/" << name; 01262 TH1F* gainChannel = new TH1F(name.c_str(),title.str().c_str(), 01263 300,0,300); 01264 gainChannel->GetXaxis()->SetTitle("Channel"); 01265 gainChannel->GetYaxis()->SetTitle("gain in ADC bin"); 01266 // DEBUG 01267 info() << "name= " << gainChannel->GetName() 01268 << " at path = \"" << path.str() << "\"" << endreq; 01269 if( m_statsSvc->put(path.str(),gainChannel).isFailure() ) { 01270 error() << "prepareStats(): Could not register " << path << endreq; 01271 delete gainChannel; 01272 return StatusCode::FAILURE; 01273 } 01274 } 01275 01276 // gain vs channel with cuts 01277 { 01278 std::ostringstream title, path; 01279 std::string name = "gainChannelWithCuts"; 01280 title << "gainChannel for readouts in " << detector.detName(); 01281 path << this->getPath(detector) << "/" << name; 01282 TH1F* gainChannelWithCuts = new TH1F(name.c_str(),title.str().c_str(), 01283 300,0,300); 01284 gainChannelWithCuts->GetXaxis()->SetTitle("channel"); 01285 gainChannelWithCuts->GetYaxis()->SetTitle("gain in ADC bin"); 01286 // DEBUG 01287 info() << "name= " << gainChannelWithCuts->GetName() 01288 << " at path = \"" << path.str() << "\"" << endreq; 01289 if( m_statsSvc->put(path.str(),gainChannelWithCuts).isFailure() ) { 01290 error() << "prepareStats(): Could not register " << path << endreq; 01291 delete gainChannelWithCuts; 01292 return StatusCode::FAILURE; 01293 } 01294 } 01295 01296 // tdc offset vs channel spectrum 01297 { 01298 std::ostringstream title, path; 01299 std::string name = "tdcOffsetChannel"; 01300 title << "tdcOffsetChannel for readouts in " << detector.detName(); 01301 path << this->getPath(detector) << "/" << name; 01302 TH1F* tdcOffsetChannel = new TH1F(name.c_str(),title.str().c_str(), 01303 300,0,300); 01304 tdcOffsetChannel->GetXaxis()->SetTitle("Channel"); 01305 tdcOffsetChannel->GetYaxis()->SetTitle("tdc offset in TDC bin"); 01306 // DEBUG 01307 info() << "name= " << tdcOffsetChannel->GetName() 01308 << " at path = \"" << path.str() << "\"" << endreq; 01309 if( m_statsSvc->put(path.str(),tdcOffsetChannel).isFailure() ) { 01310 error() << "prepareStats(): Could not register " << path << endreq; 01311 delete tdcOffsetChannel; 01312 return StatusCode::FAILURE; 01313 } 01314 } 01315 01316 // tdcOffset vs channel with cuts 01317 { 01318 std::ostringstream title, path; 01319 std::string name = "tdcOffsetChannelWithCuts"; 01320 title << "tdcOffsetChannel for readouts in " << detector.detName(); 01321 path << this->getPath(detector) << "/" << name; 01322 TH1F* tdcOffsetChannelWithCuts = new TH1F(name.c_str(),title.str().c_str(), 01323 300,0,300); 01324 tdcOffsetChannelWithCuts->GetXaxis()->SetTitle("channel"); 01325 tdcOffsetChannelWithCuts->GetYaxis()->SetTitle("tdc offset in TDC bin"); 01326 // DEBUG 01327 info() << "name= " << tdcOffsetChannelWithCuts->GetName() 01328 << " at path = \"" << path.str() << "\"" << endreq; 01329 if( m_statsSvc->put(path.str(),tdcOffsetChannelWithCuts).isFailure() ) { 01330 error() << "prepareStats(): Could not register " << path << endreq; 01331 delete tdcOffsetChannelWithCuts; 01332 return StatusCode::FAILURE; 01333 } 01334 } 01335 01336 01337 // occupancy vs channel spectrum 01338 { 01339 std::ostringstream title, path; 01340 std::string name = "occupancyChannel"; 01341 title << "occupancyChannel for readouts in " << detector.detName(); 01342 path << this->getPath(detector) << "/" << name; 01343 TH1F* occupancyChannel = new TH1F(name.c_str(),title.str().c_str(), 01344 300,0,300); 01345 occupancyChannel->GetXaxis()->SetTitle("Channel"); 01346 occupancyChannel->GetYaxis()->SetTitle("occupancy"); 01347 // DEBUG 01348 info() << "name= " << occupancyChannel->GetName() 01349 << " at path = \"" << path.str() << "\"" << endreq; 01350 if( m_statsSvc->put(path.str(),occupancyChannel).isFailure() ) { 01351 error() << "prepareStats(): Could not register " << path << endreq; 01352 delete occupancyChannel; 01353 return StatusCode::FAILURE; 01354 } 01355 } 01356 01357 // occupancy vs channel with cuts 01358 { 01359 std::ostringstream title, path; 01360 std::string name = "occupancyChannelWithCuts"; 01361 title << "occupancyChannel for readouts in " << detector.detName(); 01362 path << this->getPath(detector) << "/" << name; 01363 TH1F* occupancyChannelWithCuts = new TH1F(name.c_str(),title.str().c_str(), 01364 300,0,300); 01365 occupancyChannelWithCuts->GetXaxis()->SetTitle("channel"); 01366 occupancyChannelWithCuts->GetYaxis()->SetTitle("occupancy"); 01367 // DEBUG 01368 info() << "name= " << occupancyChannelWithCuts->GetName() 01369 << " at path = \"" << path.str() << "\"" << endreq; 01370 if( m_statsSvc->put(path.str(),occupancyChannelWithCuts).isFailure() ) { 01371 error() << "prepareStats(): Could not register " << path << endreq; 01372 delete occupancyChannelWithCuts; 01373 return StatusCode::FAILURE; 01374 } 01375 } 01376 01377 01378 // TDC Median spectrum 01379 { 01380 std::ostringstream title, path; 01381 std::string name = "tdcMedian"; 01382 title << "Median TDC for readouts in " << detector.detName(); 01383 path << this->getPath(detector) << "/" << name; 01384 TH1F* tdcMedian = new TH1F(name.c_str(),title.str().c_str(), 01385 800,0,800); 01386 tdcMedian->GetXaxis()->SetTitle("TDC value"); 01387 tdcMedian->GetYaxis()->SetTitle("Entries"); 01388 if( m_statsSvc->put(path.str(),tdcMedian).isFailure() ) { 01389 error() << "prepareStats(): Could not register " << path << endreq; 01390 delete tdcMedian; 01391 return StatusCode::FAILURE; 01392 } 01393 } 01394 01395 // TDC Median spectrum with cuts 01396 { 01397 std::ostringstream title, path; 01398 std::string name = "tdcMedianWithCuts"; 01399 title << "Median TDC WithCuts for readouts in " << detector.detName(); 01400 path << this->getPath(detector) << "/" << name; 01401 TH1F* tdcMedianWithCuts = new TH1F(name.c_str(),title.str().c_str(), 01402 800,0,800); 01403 tdcMedianWithCuts->GetXaxis()->SetTitle("TDC value"); 01404 tdcMedianWithCuts->GetYaxis()->SetTitle("Entries"); 01405 if( m_statsSvc->put(path.str(),tdcMedianWithCuts).isFailure() ) { 01406 error() << "prepareStats(): Could not register " << path << endreq; 01407 delete tdcMedianWithCuts; 01408 return StatusCode::FAILURE; 01409 } 01410 } 01411 01412 // Mean TDC Offset by AD ring 01413 { 01414 std::ostringstream title, path; 01415 std::string name = "meanTdcOffset"; 01416 title << "Mean TDC offset by ring in " << detector.detName(); 01417 path << this->getPath(detector) << "/" << name; 01418 TH1F* meanTdcOffset = new TH1F(name.c_str(),title.str().c_str(), 01419 8,1,9); 01420 meanTdcOffset->GetXaxis()->SetTitle("AD Ring"); 01421 meanTdcOffset->GetYaxis()->SetTitle("Mean TDC Offset"); 01422 if( m_statsSvc->put(path.str(),meanTdcOffset).isFailure() ) { 01423 error() << "prepareStats(): Could not register " << path << endreq; 01424 delete meanTdcOffset; 01425 return StatusCode::FAILURE; 01426 } 01427 } 01428 01429 // Mean TDC Offset WithCuts by AD ring 01430 { 01431 std::ostringstream title, path; 01432 std::string name = "meanTdcOffsetWithCuts"; 01433 title << "Mean TDC offset by ring in " << detector.detName(); 01434 path << this->getPath(detector) << "/" << name; 01435 TH1F* meanTdcOffsetWithCuts = new TH1F(name.c_str(),title.str().c_str(), 01436 8,1,9); 01437 meanTdcOffsetWithCuts->GetXaxis()->SetTitle("AD Ring"); 01438 meanTdcOffsetWithCuts->GetYaxis()->SetTitle("Mean TDC Offset With Cuts"); 01439 if( m_statsSvc->put(path.str(),meanTdcOffsetWithCuts).isFailure() ) { 01440 error() << "prepareStats(): Could not register " << path << endreq; 01441 delete meanTdcOffsetWithCuts; 01442 return StatusCode::FAILURE; 01443 } 01444 } 01445 01446 // Mean Occupancy by AD ring 01447 { 01448 std::ostringstream title, path; 01449 std::string name = "meanOccupancy"; 01450 title << "Mean occupancy by ring in " << detector.detName(); 01451 path << this->getPath(detector) << "/" << name; 01452 TH1F* meanOccupancy = new TH1F(name.c_str(),title.str().c_str(), 01453 8,1,9); 01454 meanOccupancy->GetXaxis()->SetTitle("AD Ring"); 01455 meanOccupancy->GetYaxis()->SetTitle("Mean Occupancy"); 01456 if( m_statsSvc->put(path.str(),meanOccupancy).isFailure() ) { 01457 error() << "prepareStats(): Could not register " << path << endreq; 01458 delete meanOccupancy; 01459 return StatusCode::FAILURE; 01460 } 01461 } 01462 01463 // Mean Adc by AD ring 01464 { 01465 std::ostringstream title, path; 01466 std::string name = "meanAdc"; 01467 title << "Mean Adc by ring in " << detector.detName(); 01468 path << this->getPath(detector) << "/" << name; 01469 TH1F* meanAdc = new TH1F(name.c_str(),title.str().c_str(), 01470 8,1,9); 01471 meanAdc->GetXaxis()->SetTitle("AD Ring"); 01472 meanAdc->GetYaxis()->SetTitle("Mean Adc"); 01473 if( m_statsSvc->put(path.str(),meanAdc).isFailure() ) { 01474 error() << "prepareStats(): Could not register " << path << endreq; 01475 delete meanAdc; 01476 return StatusCode::FAILURE; 01477 } 01478 } 01479 01480 // Mean Occupancy With Cuts by AD ring 01481 { 01482 std::ostringstream title, path; 01483 std::string name = "meanOccupancyWithCuts"; 01484 title << "Mean occupancy With Cuts by ring in " << detector.detName(); 01485 path << this->getPath(detector) << "/" << name; 01486 TH1F* meanOccupancyWithCuts = new TH1F(name.c_str(),title.str().c_str(), 01487 8,1,9); 01488 meanOccupancyWithCuts->GetXaxis()->SetTitle("AD Ring"); 01489 meanOccupancyWithCuts->GetYaxis()->SetTitle("Mean Occupancy With Cuts"); 01490 if( m_statsSvc->put(path.str(),meanOccupancyWithCuts).isFailure() ) { 01491 error() << "prepareStats(): Could not register " << path << endreq; 01492 delete meanOccupancyWithCuts; 01493 return StatusCode::FAILURE; 01494 } 01495 } 01496 01497 // Mean Adc With Cuts by AD ring 01498 { 01499 std::ostringstream title, path; 01500 std::string name = "meanAdcWithCuts"; 01501 title << "Mean Adc With Cuts by ring in " << detector.detName(); 01502 path << this->getPath(detector) << "/" << name; 01503 TH1F* meanAdcWithCuts = new TH1F(name.c_str(),title.str().c_str(), 01504 8,1,9); 01505 meanAdcWithCuts->GetXaxis()->SetTitle("AD Ring"); 01506 meanAdcWithCuts->GetYaxis()->SetTitle("Mean Adc With Cuts"); 01507 if( m_statsSvc->put(path.str(),meanAdcWithCuts).isFailure() ) { 01508 error() << "prepareStats(): Could not register " << path << endreq; 01509 delete meanAdcWithCuts; 01510 return StatusCode::FAILURE; 01511 } 01512 } 01513 01514 01515 // Make sure summary histograms for each ring exist in detector statistics 01516 for(int ring = 0; ring <= 8; ring++){ 01517 // Occupancy by ring 01518 { 01519 std::ostringstream title, path; 01520 std::string name = "occupancy"; 01521 title << "Occupancy for PMTs in " << detector.detName() 01522 << " ring " << ring; 01523 path << this->getPath(detector, ring) << "/" << name; 01524 TH1F* occupancy = new TH1F(name.c_str(),title.str().c_str(), 01525 24,1,25); 01526 occupancy->GetXaxis()->SetTitle("Column"); 01527 occupancy->GetYaxis()->SetTitle("Occupancy"); 01528 if( m_statsSvc->put(path.str(),occupancy).isFailure() ) { 01529 error() << "prepareStats(): Could not register " << path << endreq; 01530 delete occupancy; 01531 return StatusCode::FAILURE; 01532 } 01533 } 01534 // TDC offset by ring 01535 { 01536 std::ostringstream title, path; 01537 std::string name = "tdcOffset"; 01538 title << "Time Offset for PMTs in " << detector.detName() 01539 << " ring " << ring; 01540 path << this->getPath(detector, ring) << "/" << name; 01541 TH1F* tdcOffset = new TH1F(name.c_str(),title.str().c_str(), 01542 24,1,25); 01543 tdcOffset->GetXaxis()->SetTitle("Column"); 01544 tdcOffset->GetYaxis()->SetTitle("Time Offset [tdc]"); 01545 if( m_statsSvc->put(path.str(),tdcOffset).isFailure() ) { 01546 error() << "prepareStats(): Could not register " << path << endreq; 01547 delete tdcOffset; 01548 return StatusCode::FAILURE; 01549 } 01550 } 01551 // TDC Raw offset by ring 01552 { 01553 std::ostringstream title, path; 01554 std::string name = "tdcRawOffset"; 01555 title << "Time Raw Offset for PMTs in " << detector.detName() 01556 << " ring " << ring; 01557 path << this->getPath(detector, ring) << "/" << name; 01558 TH1F* tdcRawOffset = new TH1F(name.c_str(),title.str().c_str(), 01559 24,1,25); 01560 tdcRawOffset->GetXaxis()->SetTitle("Column"); 01561 tdcRawOffset->GetYaxis()->SetTitle("Time Offset [tdc]"); 01562 if( m_statsSvc->put(path.str(),tdcRawOffset).isFailure() ) { 01563 error() << "prepareStats(): Could not register " << path << endreq; 01564 delete tdcRawOffset; 01565 return StatusCode::FAILURE; 01566 } 01567 } 01568 01569 // ADC Median by ring 01570 { 01571 std::ostringstream title, path; 01572 std::string name = "adcMedian"; 01573 title << "ADC Median for PMTs in " << detector.detName() 01574 << " ring " << ring; 01575 path << this->getPath(detector, ring) << "/" << name; 01576 TH1F* adcMedian = new TH1F(name.c_str(),title.str().c_str(), 01577 24,1,25); 01578 adcMedian->GetXaxis()->SetTitle("Column"); 01579 adcMedian->GetYaxis()->SetTitle("Median ADC"); 01580 if( m_statsSvc->put(path.str(),adcMedian).isFailure() ) { 01581 error() << "prepareStats(): Could not register " << path << endreq; 01582 delete adcMedian; 01583 return StatusCode::FAILURE; 01584 } 01585 } 01586 // ADC Sigma by ring 01587 { 01588 std::ostringstream title, path; 01589 std::string name = "adcSigma"; 01590 title << "ADC Sigma for PMTs in " << detector.detName() 01591 << " ring " << ring; 01592 path << this->getPath(detector, ring) << "/" << name; 01593 TH1F* adcSigma = new TH1F(name.c_str(),title.str().c_str(), 01594 24,1,25); 01595 adcSigma->GetXaxis()->SetTitle("Column"); 01596 adcSigma->GetYaxis()->SetTitle("ADC Sigma"); 01597 if( m_statsSvc->put(path.str(),adcSigma).isFailure() ) { 01598 error() << "prepareStats(): Could not register " << path << endreq; 01599 delete adcSigma; 01600 return StatusCode::FAILURE; 01601 } 01602 } 01603 } 01604 01605 // Make sure summary histograms for each ring exist in detector statistics with cuts 01606 for(int ring = 0; ring <= 8; ring++){ 01607 // Occupancy by ring With Cuts 01608 { 01609 std::ostringstream title, path; 01610 std::string name = "occupancyWithCuts"; 01611 title << "Occupancy With Cuts for PMTs in " << detector.detName() 01612 << " ring " << ring; 01613 path << this->getPath(detector, ring) << "/" << name; 01614 TH1F* occupancyWithCuts = new TH1F(name.c_str(),title.str().c_str(), 01615 24,1,25); 01616 occupancyWithCuts->GetXaxis()->SetTitle("Column"); 01617 occupancyWithCuts->GetYaxis()->SetTitle("Occupancy"); 01618 if( m_statsSvc->put(path.str(),occupancyWithCuts).isFailure() ) { 01619 error() << "prepareStats(): Could not register " << path << endreq; 01620 delete occupancyWithCuts; 01621 return StatusCode::FAILURE; 01622 } 01623 } 01624 // TDC offset by ring With Cuts 01625 { 01626 std::ostringstream title, path; 01627 std::string name = "tdcOffsetWithCuts"; 01628 title << "Time Offset With Cuts for PMTs in " << detector.detName() 01629 << " ring " << ring; 01630 path << this->getPath(detector, ring) << "/" << name; 01631 TH1F* tdcOffsetWithCuts = new TH1F(name.c_str(),title.str().c_str(), 01632 24,1,25); 01633 tdcOffsetWithCuts->GetXaxis()->SetTitle("Column"); 01634 tdcOffsetWithCuts->GetYaxis()->SetTitle("Time Offset With Cuts [tdc]"); 01635 if( m_statsSvc->put(path.str(),tdcOffsetWithCuts).isFailure() ) { 01636 error() << "prepareStats(): Could not register " << path << endreq; 01637 delete tdcOffsetWithCuts; 01638 return StatusCode::FAILURE; 01639 } 01640 } 01641 // ADC Median by ring With Cuts 01642 { 01643 std::ostringstream title, path; 01644 std::string name = "adcMedianWithCuts"; 01645 title << "ADC Median With Cuts for PMTs in " << detector.detName() 01646 << " ring " << ring; 01647 path << this->getPath(detector, ring) << "/" << name; 01648 TH1F* adcMedianWithCuts = new TH1F(name.c_str(),title.str().c_str(), 01649 24,1,25); 01650 adcMedianWithCuts->GetXaxis()->SetTitle("Column"); 01651 adcMedianWithCuts->GetYaxis()->SetTitle("Median ADC With Cuts"); 01652 if( m_statsSvc->put(path.str(),adcMedianWithCuts).isFailure() ) { 01653 error() << "prepareStats(): Could not register " << path << endreq; 01654 delete adcMedianWithCuts; 01655 return StatusCode::FAILURE; 01656 } 01657 } 01658 // ADC Sigma by ring With Cuts 01659 { 01660 std::ostringstream title, path; 01661 std::string name = "adcSigmaWithCuts"; 01662 title << "ADC Sigma With Cuts for PMTs in " << detector.detName() 01663 << " ring " << ring; 01664 path << this->getPath(detector, ring) << "/" << name; 01665 TH1F* adcSigmaWithCuts = new TH1F(name.c_str(),title.str().c_str(), 01666 24,1,25); 01667 adcSigmaWithCuts->GetXaxis()->SetTitle("Column"); 01668 adcSigmaWithCuts->GetYaxis()->SetTitle("ADC Sigma With Cuts"); 01669 if( m_statsSvc->put(path.str(),adcSigmaWithCuts).isFailure() ) { 01670 error() << "prepareStats(): Could not register " << path << endreq; 01671 delete adcSigmaWithCuts; 01672 return StatusCode::FAILURE; 01673 } 01674 } 01675 } 01676 01677 01678 ServiceMode svcMode(context, 0); 01679 01680 // Get list of all FEE channels 01681 const std::vector<DayaBay::FeeChannelId>& channelList 01682 = m_cableSvc->feeChannelIds( svcMode ); 01683 std::vector<DayaBay::FeeChannelId>::const_iterator chanIter, 01684 chanEnd = channelList.end(); 01685 // Initialize statistics for each channel 01686 for(chanIter = channelList.begin(); chanIter != chanEnd; chanIter++){ 01687 DayaBay::FeeChannelId chanId = *chanIter; 01688 01689 // Board Number 01690 { 01691 std::string name = "board"; 01692 std::ostringstream path; 01693 path << this->getPath(chanId) << "/" << name; 01694 TParameter<int>* par = new TParameter<int>(name.c_str(), chanId.board()); 01695 if( m_statsSvc->put(path.str(),par).isFailure() ) { 01696 error() << "prepareStats(): Could not register " << name 01697 << " at " << path << endreq; 01698 delete par; 01699 par = 0; 01700 return StatusCode::FAILURE; 01701 } 01702 } 01703 01704 // Connector Number 01705 { 01706 std::string name = "connector"; 01707 std::ostringstream path; 01708 path << this->getPath(chanId) << "/" << name; 01709 TParameter<int>* par = new TParameter<int>(name.c_str(), 01710 chanId.connector()); 01711 if( m_statsSvc->put(path.str(),par).isFailure() ) { 01712 error() << "prepareStats(): Could not register " << name 01713 << " at " << path << endreq; 01714 delete par; 01715 par = 0; 01716 return StatusCode::FAILURE; 01717 } 01718 } 01719 01720 // Number of Hits 01721 { 01722 std::string name = "nHits"; 01723 std::ostringstream path; 01724 path << this->getPath(chanId) << "/" << name; 01725 TParameter<int>* par = new TParameter<int>(name.c_str(), 0); 01726 if( m_statsSvc->put(path.str(),par).isFailure() ) { 01727 error() << "prepareStats(): Could not register " << name 01728 << " at " << path << endreq; 01729 delete par; 01730 par = 0; 01731 return StatusCode::FAILURE; 01732 } 01733 } 01734 01735 // Number of Hits with cuts 01736 { 01737 std::string name = "nHitsWithCuts"; 01738 std::ostringstream path; 01739 path << this->getPath(chanId) << "/" << name; 01740 TParameter<int>* parWithCuts = new TParameter<int>(name.c_str(), 0); 01741 if( m_statsSvc->put(path.str(),parWithCuts).isFailure() ) { 01742 error() << "prepareStats(): Could not register " << name 01743 << " at " << path << endreq; 01744 delete parWithCuts; 01745 parWithCuts = 0; 01746 return StatusCode::FAILURE; 01747 } 01748 } 01749 01750 01751 // Raw TDC spectrum 01752 { 01753 std::string name = "tdcRaw"; 01754 std::ostringstream title, path; 01755 title << "Raw TDC values for " << detector.detName() 01756 << " channel " << chanId.board() << "_" 01757 << chanId.connector(); 01758 path << this->getPath(chanId) << "/" << name; 01759 TH1F* tdcRaw = new TH1F(name.c_str(),title.str().c_str(), 01760 300,0,300); 01761 tdcRaw->GetXaxis()->SetTitle("TDC value"); 01762 tdcRaw->GetYaxis()->SetTitle("Entries"); 01763 if( m_statsSvc->put(path.str(),tdcRaw).isFailure() ) { 01764 error() << "prepareStats(): Could not register " << path << endreq; 01765 delete tdcRaw; 01766 return StatusCode::FAILURE; 01767 } 01768 } 01769 01770 // TDC spectrum with cuts 01771 { 01772 std::string name = "tdcWithCuts"; 01773 std::ostringstream title, path; 01774 title << "TDC values With Cuts for " << detector.detName() 01775 << " channel " << chanId.board() << "_" 01776 << chanId.connector(); 01777 path << this->getPath(chanId) << "/" << name; 01778 TH1F* tdcWithCuts = new TH1F(name.c_str(),title.str().c_str(), 01779 800,0,800); 01780 tdcWithCuts->GetXaxis()->SetTitle("TDC value"); 01781 tdcWithCuts->GetYaxis()->SetTitle("Entries With Cuts"); 01782 if( m_statsSvc->put(path.str(),tdcWithCuts).isFailure() ) { 01783 error() << "prepareStats(): Could not register " << path << endreq; 01784 delete tdcWithCuts; 01785 return StatusCode::FAILURE; 01786 } 01787 } 01788 01789 01790 // TDC spectrum, median corrected 01791 { 01792 std::string name = "tdcByMedian"; 01793 std::ostringstream title, path; 01794 title << "Median-corrected TDC values for " << detector.detName() 01795 << " channel " << chanId.board() << "_" 01796 << chanId.connector(); 01797 path << this->getPath(chanId) << "/" << name; 01798 TH1F* tdcByMedian = new TH1F(name.c_str(),title.str().c_str(), 01799 600,-300,300); 01800 tdcByMedian->GetXaxis()->SetTitle("TDC value"); 01801 tdcByMedian->GetYaxis()->SetTitle("Entries"); 01802 if( m_statsSvc->put(path.str(),tdcByMedian).isFailure() ) { 01803 error() << "prepareStats(): Could not register " << path << endreq; 01804 delete tdcByMedian; 01805 return StatusCode::FAILURE; 01806 } 01807 } 01808 01809 // TDC spectrum, median corrected with cuts 01810 { 01811 std::string name = "tdcByMedianWithCuts"; 01812 std::ostringstream title, path; 01813 title << "Median-corrected TDC values With Cuts for " << detector.detName() 01814 << " channel " << chanId.board() << "_" 01815 << chanId.connector(); 01816 path << this->getPath(chanId) << "/" << name; 01817 TH1F* tdcByMedianWithCuts = new TH1F(name.c_str(),title.str().c_str(), 01818 600,-300,300); 01819 tdcByMedianWithCuts->GetXaxis()->SetTitle("TDC value"); 01820 tdcByMedianWithCuts->GetYaxis()->SetTitle("Entries With Cuts"); 01821 if( m_statsSvc->put(path.str(),tdcByMedianWithCuts).isFailure() ) { 01822 error() << "prepareStats(): Could not register " << path << endreq; 01823 delete tdcByMedianWithCuts; 01824 return StatusCode::FAILURE; 01825 } 01826 } 01827 01828 01829 // Raw ADC spectrum 01830 { 01831 std::ostringstream title, path; 01832 std::string name = "adcRaw"; 01833 title << "ADC values for " << detector.detName() 01834 << " channel " << chanId.board() << "_" 01835 << chanId.connector(); 01836 path << this->getPath(chanId) << "/" << name; 01837 TH1F* adcRaw = new TH1F(name.c_str(),title.str().c_str(), 01838 4096,0,4096); 01839 adcRaw->GetXaxis()->SetTitle("ADC value"); 01840 adcRaw->GetYaxis()->SetTitle("Entries"); 01841 if( m_statsSvc->put(path.str(),adcRaw).isFailure() ) { 01842 error() << "prepareStats(): Could not register " << path << endreq; 01843 delete adcRaw; 01844 return StatusCode::FAILURE; 01845 } 01846 } 01847 01848 // ADC spectrum with cuts 01849 { 01850 std::ostringstream title, path; 01851 std::string name = "adcRawWithCuts"; 01852 title << "ADC values for " << detector.detName() 01853 << " channel " << chanId.board() << "_" 01854 << chanId.connector(); 01855 path << this->getPath(chanId) << "/" << name; 01856 TH1F* adcRawWithCuts = new TH1F(name.c_str(),title.str().c_str(), 01857 4096,0,4096); 01858 adcRawWithCuts->GetXaxis()->SetTitle("ADC value"); 01859 adcRawWithCuts->GetYaxis()->SetTitle("Entries With Cuts"); 01860 if( m_statsSvc->put(path.str(),adcRawWithCuts).isFailure() ) { 01861 error() << "prepareStats(): Could not register " << path << endreq; 01862 delete adcRawWithCuts; 01863 return StatusCode::FAILURE; 01864 } 01865 } 01866 01867 01868 // ADC spectrum, baseline subtracted 01869 { 01870 std::ostringstream title, path; 01871 std::string name = "adc"; 01872 title << "Baseline-subtracted ADC values for " << detector.detName() 01873 << " channel " << chanId.board() << "_" 01874 << chanId.connector(); 01875 path << this->getPath(chanId) << "/" << name; 01876 TH1F* adc = new TH1F(name.c_str(),title.str().c_str(), 01877 5096,-1000,4096); 01878 adc->GetXaxis()->SetTitle("ADC value"); 01879 adc->GetYaxis()->SetTitle("Entries"); 01880 if( m_statsSvc->put(path.str(),adc).isFailure() ) { 01881 error() << "prepareStats(): Could not register " << path << endreq; 01882 delete adc; 01883 return StatusCode::FAILURE; 01884 } 01885 } 01886 01887 // ADC spectrum, baseline subtracted with cuts 01888 { 01889 std::ostringstream title, path; 01890 std::string name = "adcWithCuts"; 01891 title << "Baseline-subtracted ADC values With Cuts for " << detector.detName() 01892 << " channel " << chanId.board() << "_" 01893 << chanId.connector(); 01894 path << this->getPath(chanId) << "/" << name; 01895 TH1F* adcWithCuts = new TH1F(name.c_str(),title.str().c_str(), 01896 5096,-1000,4096); 01897 adcWithCuts->GetXaxis()->SetTitle("ADC value"); 01898 adcWithCuts->GetYaxis()->SetTitle("Entries"); 01899 if( m_statsSvc->put(path.str(),adcWithCuts).isFailure() ) { 01900 error() << "prepareStats(): Could not register " << path << endreq; 01901 delete adcWithCuts; 01902 return StatusCode::FAILURE; 01903 } 01904 } 01905 01906 01907 // ADC spectrum by Clock Cycle 01908 { 01909 std::ostringstream title, path; 01910 std::string name = "adcByClock"; 01911 title << "ADC values by clock cycle for " << detector.detName() 01912 << " channel " << chanId.board() << "_" 01913 << chanId.connector(); 01914 path << this->getPath(chanId) << "/" << name; 01915 TH1F* adcByClock = new TH1F(name.c_str(),title.str().c_str(), 01916 20,0,20); 01917 adcByClock->GetXaxis()->SetTitle("ADC Clock Cycle"); 01918 adcByClock->GetYaxis()->SetTitle("ADC"); 01919 if( m_statsSvc->put(path.str(),adcByClock).isFailure() ) { 01920 error() << "prepareStats(): Could not register " << path << endreq; 01921 delete adcByClock; 01922 return StatusCode::FAILURE; 01923 } 01924 } 01925 // ADC spectrum by Clock Cycle with cuts 01926 { 01927 std::ostringstream title, path; 01928 std::string name = "adcByClockWithCuts"; 01929 title << "ADC values by clock cycle for " << detector.detName() 01930 << " channel " << chanId.board() << "_" 01931 << chanId.connector(); 01932 path << this->getPath(chanId) << "/" << name; 01933 TH1F* adcByClockWithCuts = new TH1F(name.c_str(),title.str().c_str(), 01934 20,0,20); 01935 adcByClockWithCuts->GetXaxis()->SetTitle("ADC Clock Cycle"); 01936 adcByClockWithCuts->GetYaxis()->SetTitle("ADC"); 01937 if( m_statsSvc->put(path.str(),adcByClockWithCuts).isFailure() ) { 01938 error() << "prepareStats(): Could not register " << path << endreq; 01939 delete adcByClockWithCuts; 01940 return StatusCode::FAILURE; 01941 } 01942 } 01943 01944 01945 } 01946 01947 m_processedDetectors.push_back(detector); 01948 01949 return StatusCode::SUCCESS; 01950 }
| std::string PmtCalibLeadingEdgeWithCuts::getPath | ( | const DayaBay::Detector & | detector | ) | [private] |
Definition at line 1952 of file PmtCalibLeadingEdgeWithCuts.cc.
01952 { 01953 return m_filepath + "/" + detector.detName(); 01954 }
| std::string PmtCalibLeadingEdgeWithCuts::getPath | ( | const DayaBay::Detector & | detector, | |
| int | ring | |||
| ) | [private] |
Definition at line 1956 of file PmtCalibLeadingEdgeWithCuts.cc.
01957 { 01958 std::ostringstream path; 01959 path << m_filepath << "/" << detector.detName() << "/ring_" << ring; 01960 return path.str(); 01961 }
| std::string PmtCalibLeadingEdgeWithCuts::getPath | ( | const DayaBay::FeeChannelId & | channelId | ) | [private] |
Definition at line 1963 of file PmtCalibLeadingEdgeWithCuts.cc.
01963 { 01964 std::ostringstream path; 01965 path << m_filepath << "/" << chanId.detName() 01966 << "/chan_" << chanId.board() << "_" << chanId.connector(); 01967 return path.str(); 01968 }
| StatusCode PmtCalibLeadingEdgeWithCuts::calibrateRaw | ( | ) | [private, virtual] |
used for firstly calibration of PMT before calibration with cuts calibrateRaw() is called after processing many readouts.
This method is responsible for calculating the calibration parameters. wangzhm@ihep.ac.cn 2009/11/28
Definition at line 391 of file PmtCalibLeadingEdgeWithCuts.cc.
00391 { 00392 // Loop over detectors in summary data, and calculate parameters 00393 00394 std::vector<DayaBay::Detector>::iterator detIter, 00395 detEnd = m_processedDetectors.end(); 00396 for(detIter = m_processedDetectors.begin(); detIter != detEnd; detIter++){ 00397 DayaBay::Detector detector = *detIter; 00398 00399 // Retrieve the data description and histograms 00400 TObject* nReadoutsObj = 00401 m_statsSvc->get(this->getPath(detector) + "/nReadouts"); 00402 TParameter<int>* nReadoutsPar= dynamic_cast<TParameter<int>*>(nReadoutsObj); 00403 if(!nReadoutsPar) return StatusCode::FAILURE; 00404 00405 TH1F* meanOccupancyH = m_statsSvc->getTH1F( this->getPath(detector) 00406 +"/meanOccupancy"); 00407 if(!meanOccupancyH) return StatusCode::FAILURE; 00408 00409 TH1F* meanAdcH = m_statsSvc->getTH1F( this->getPath(detector) 00410 +"/meanAdc"); 00411 if(!meanAdcH) return StatusCode::FAILURE; 00412 00413 TH1F* gainChannelH = m_statsSvc->getTH1F( this->getPath(detector) 00414 +"/gainChannel"); 00415 if(!gainChannelH) return StatusCode::FAILURE; 00416 TH1F* tdcOffsetChannelH = m_statsSvc->getTH1F( this->getPath(detector) 00417 +"/tdcOffsetChannel"); 00418 if(!tdcOffsetChannelH) return StatusCode::FAILURE; 00419 TH1F* occupancyChannelH = m_statsSvc->getTH1F( this->getPath(detector) 00420 +"/occupancyChannel"); 00421 if(!occupancyChannelH) return StatusCode::FAILURE; 00422 00423 00424 TH1F* meanTdcOffsetH = m_statsSvc->getTH1F( this->getPath(detector) 00425 +"/meanTdcOffset"); 00426 if(!meanTdcOffsetH) return StatusCode::FAILURE; 00427 00428 TObject* simFlagObj = m_statsSvc->get(this->getPath(detector) +"/simFlag"); 00429 TParameter<int>* simFlagPar = dynamic_cast<TParameter<int>*>(simFlagObj); 00430 if(!simFlagPar) return StatusCode::FAILURE; 00431 00432 TObject* timeSecObj = m_statsSvc->get(this->getPath(detector) +"/timeSec"); 00433 TParameter<int>* timeSecPar = dynamic_cast<TParameter<int>*>(timeSecObj); 00434 if(!timeSecPar) return StatusCode::FAILURE; 00435 00436 TObject* timeNanoSecObj = 00437 m_statsSvc->get(this->getPath(detector) +"/timeNanoSec"); 00438 TParameter<int>* timeNanoSecPar = 00439 dynamic_cast<TParameter<int>*>(timeNanoSecObj); 00440 if(!timeNanoSecPar) return StatusCode::FAILURE; 00441 00442 Site::Site_t site = detector.site(); 00443 DetectorId::DetectorId_t detId = detector.detectorId(); 00444 SimFlag::SimFlag_t simFlag = (SimFlag::SimFlag_t)(simFlagPar->GetVal()); 00445 time_t timeSec = (time_t)(timeSecPar->GetVal()); 00446 int timeNanoSec = timeNanoSecPar->GetVal(); 00447 int nReadouts = nReadoutsPar->GetVal(); 00448 00449 // Context, ServiceMode 00450 Context context(site,simFlag,TimeStamp(timeSec,timeNanoSec),detId); 00451 ServiceMode svcMode(context, 0); 00452 00453 // Prepare data for ring-to-ring comparison 00454 int nRings = 8; 00455 std::map<int,double> meanOccRing; 00456 std::map<int,double> meanOccWeight; 00457 std::map<int,double> meanTdcRing; 00458 std::map<int,double> meanTdcWeight; 00459 for(int ring = 1; ring <= nRings; ring++){ 00460 meanOccRing[ring] = 0.0; 00461 meanOccWeight[ring] = 0.0; 00462 meanTdcRing[ring] = 0.0; 00463 meanTdcWeight[ring] = 0.0; 00464 } 00465 double minValidOcc = 0.01; // Cut on PMTs with bad/low occupancy 00466 std::vector<DayaBay::FeeChannelId> badChannels; 00467 00468 // Get list of all FEE channels 00469 const std::vector<DayaBay::FeeChannelId>& channelList 00470 = m_cableSvc->feeChannelIds( svcMode ); 00471 std::vector<DayaBay::FeeChannelId>::const_iterator chanIter, 00472 chanEnd = channelList.end(); 00473 // Initialize statistics for each channel 00474 for(chanIter = channelList.begin(); chanIter != chanEnd; chanIter++){ 00475 DayaBay::FeeChannelId chanId = *chanIter; 00476 00477 // Analyze the channel data 00478 // Determine ring and column 00479 DayaBay::AdPmtSensor pmtId = 00480 m_cableSvc->adPmtSensor(chanId, svcMode); 00481 int ring = pmtId.ring(); 00482 int column = pmtId.column(); 00483 int channel_number=ring*100+column; 00484 00485 // Channel status 00486 bool occupancyOk = true; 00487 bool gainOk = true; 00488 00489 // Get Parameters and Histograms 00490 TObject* nHitsObj = m_statsSvc->get(this->getPath(chanId) + "/nHits"); 00491 TParameter<int>* nHitsPar = dynamic_cast<TParameter<int>*>(nHitsObj); 00492 if(!nHitsPar) return StatusCode::FAILURE; 00493 TH1F* tdcByMedianH= m_statsSvc->getTH1F( this->getPath(chanId) 00494 +"/tdcByMedian"); 00495 if(!tdcByMedianH) return StatusCode::FAILURE; 00496 TH1F* tdcRawH= m_statsSvc->getTH1F( this->getPath(chanId) 00497 +"/tdcRaw"); 00498 if(!tdcRawH) return StatusCode::FAILURE; 00499 TH1F* adcH = m_statsSvc->getTH1F( this->getPath(chanId) + "/adc"); 00500 if(!adcH) return StatusCode::FAILURE; 00501 TH1F* occupancyH = m_statsSvc->getTH1F( this->getPath(detector,ring) 00502 +"/occupancy"); 00503 if(!occupancyH) return StatusCode::FAILURE; 00504 TH1F* tdcOffsetH = m_statsSvc->getTH1F( this->getPath(detector,ring) 00505 +"/tdcOffset"); 00506 if(!tdcOffsetH) return StatusCode::FAILURE; 00507 TH1F* tdcRawOffsetH = m_statsSvc->getTH1F( this->getPath(detector,ring) 00508 +"/tdcRawOffset"); 00509 if(!tdcRawOffsetH) return StatusCode::FAILURE; 00510 TH1F* adcMedianH = m_statsSvc->getTH1F( this->getPath(detector,ring) 00511 +"/adcMedian"); 00512 if(!adcMedianH) return StatusCode::FAILURE; 00513 TH1F* adcSigmaH = m_statsSvc->getTH1F( this->getPath(detector,ring) 00514 +"/adcSigma"); 00515 if(!adcSigmaH) return StatusCode::FAILURE; 00516 00517 // Channel Occupancy 00518 double occupancy = 0; 00519 double occupancyUncert = 0; 00520 { 00521 int nHits = nHitsPar->GetVal(); 00522 if(nReadouts > 0){ 00523 occupancy = nHits / ((double) nReadouts); 00524 occupancyUncert = 00525 sqrt( occupancy*(1-occupancy)/((double) nReadouts) ); 00526 } 00527 occupancyH->SetBinContent(occupancyH->FindBin(column),occupancy); 00528 occupancyH->SetBinError(occupancyH->FindBin(column),occupancyUncert); 00529 if(occupancy < minValidOcc) occupancyOk = false; 00530 00531 occupancyChannelH->SetBinContent(ring*24+column,occupancy); 00532 occupancyChannelH->SetBinError(ring*24+column,occupancyUncert); 00533 00534 } 00535 00536 // TDC offset by Median 00537 double tdcOffset = 0; 00538 double tdcOffsetUncert = 0; 00539 { 00540 // Use a fit to the leading edge to determine the time offset 00541 int maxBin = tdcByMedianH->GetMaximumBin(); 00542 double fitMin = tdcByMedianH->GetBinCenter(maxBin - 2); 00543 double fitMax = tdcByMedianH->GetBinCenter(maxBin + 10); 00544 tdcByMedianH->Fit("gaus","","",fitMin, fitMax); 00545 TF1* fitFunc = tdcByMedianH->GetFunction("gaus"); 00546 if( fitFunc ){ 00547 tdcOffset = fitFunc->GetParameter(1); 00548 tdcOffsetUncert = fitFunc->GetParError(1); 00549 00550 if( fitFunc->GetNDF() < 3 ) 00551 // Catch bad fits 00552 tdcOffsetUncert = 0; 00553 else 00554 // Scale uncertainty by the quality of the fit 00555 tdcOffsetUncert *= sqrt(fitFunc->GetChisquare() / fitFunc->GetNDF()); 00556 00557 tdcOffsetH->SetBinContent(tdcOffsetH->FindBin(column),tdcOffset); 00558 tdcOffsetH->SetBinError(tdcOffsetH->FindBin(column),tdcOffsetUncert); 00559 00560 tdcOffsetChannelH->SetBinContent(ring*24+column,tdcOffset); 00561 tdcOffsetChannelH->SetBinError(ring*24+column,tdcOffsetUncert); 00562 00563 }else{ 00564 warning() << "PMT " << pmtId << " has no TDC data for fitting." 00565 << endreq; 00566 } 00567 00568 } 00569 00570 // TDC offset Raw 00571 tdcOffset = 0; 00572 tdcOffsetUncert = 0; 00573 { 00574 // Use a fit to the leading edge to determine the time offset 00575 int maxBin = tdcRawH->GetMaximumBin(); 00576 double fitMin = tdcRawH->GetBinCenter(maxBin - 2); 00577 double fitMax = tdcRawH->GetBinCenter(maxBin + 10); 00578 tdcRawH->Fit("gaus","","",fitMin, fitMax); 00579 TF1* fitFunc = tdcRawH->GetFunction("gaus"); 00580 if( fitFunc ){ 00581 tdcOffset = fitFunc->GetParameter(1); 00582 tdcOffsetUncert = fitFunc->GetParError(1); 00583 00584 if( fitFunc->GetNDF() < 3 ) 00585 // Catch bad fits 00586 tdcOffsetUncert = 0; 00587 else 00588 // Scale uncertainty by the quality of the fit 00589 tdcOffsetUncert *= sqrt(fitFunc->GetChisquare() / fitFunc->GetNDF()); 00590 00591 tdcRawOffsetH->SetBinContent(tdcOffsetH->FindBin(column),tdcOffset); 00592 tdcRawOffsetH->SetBinError(tdcOffsetH->FindBin(column),tdcOffsetUncert); 00593 }else{ 00594 warning() << "PMT " << pmtId << " has no TDC data for fitting." 00595 << endreq; 00596 } 00597 //set for cuts of TDC 00598 if(m_TDC_Cut_Min>=0) 00599 { 00600 m_TDC_Cut_Min_Channel[channel_number]=tdcOffset-m_TDC_Cut_Min; 00601 } 00602 else 00603 { 00604 m_TDC_Cut_Min_Channel[channel_number]=0; 00605 } 00606 if(m_TDC_Cut_Max>=0) 00607 { 00608 m_TDC_Cut_Max_Channel[channel_number]=tdcOffset+m_TDC_Cut_Max; 00609 } 00610 else 00611 { 00612 m_TDC_Cut_Max_Channel[channel_number]=10000; 00613 } 00614 00615 //info() << "Issues with fit on ring, column: " << ring << " / " << column << " m_TDC_Cut_Min_Channel/m_TDC_Cut_Max_Channel: " << m_TDC_Cut_Min_Channel[channel_number] << " / " << m_TDC_Cut_Max_Channel[channel_number] << endreq; 00616 00617 } 00618 00619 00620 00621 // ADC median and sigma 00622 double adcMean = 0; 00623 double adcMeanUncert = 0; 00624 double adcSigma = 0; 00625 double adcSigmaUncert = 0; 00626 const DayaBay::FeeCalibData* feeCalib = m_calibSvc->feeCalibData(chanId, svcMode); 00627 double adcBaseline = feeCalib->m_adcBaselineLow; 00628 { 00629 // Find ADC SPE peak, width with simple gaussian 00630 TF1 adcFit("adcFit","gaus"); 00631 double adcPeak = adcH->GetBinCenter(adcH->GetMaximumBin()); 00632 double adcSig = 5.0; 00633 adcH->Fit(&adcFit,"","",adcPeak - adcSig, adcPeak + adcSig); 00634 adcMean = adcFit.GetParameter(1); 00635 if(fabs(adcMean-adcPeak) > adcSig){ 00636 warning() << "Issues with fit on ring, column: " << ring << " / " 00637 << column << " mean/peak: " << adcMean << " / " << adcPeak 00638 << endreq; 00639 gainOk = false; 00640 } 00641 //set for cuts of ADC 00642 if(m_ADC_Cut_Min>=0) 00643 { 00644 m_ADC_Cut_Min_Channel[channel_number]=(adcMean-m_ADC_Cut_Min*(adcMean-adcBaseline)); 00645 } 00646 else 00647 { 00648 m_ADC_Cut_Min_Channel[channel_number]=0; 00649 } 00650 if(m_ADC_Cut_Max>=0) 00651 { 00652 m_ADC_Cut_Max_Channel[channel_number]=(adcMean+m_ADC_Cut_Max*(adcMean-adcBaseline)); 00653 } 00654 else 00655 { 00656 m_ADC_Cut_Max_Channel[channel_number]=4096; 00657 } 00658 00659 adcMeanUncert = adcFit.GetParError(1); 00660 adcSigma = adcFit.GetParameter(2); 00661 adcSigmaUncert = adcFit.GetParError(2); 00662 adcMedianH->SetBinContent(adcMedianH->FindBin(column),adcMean); 00663 adcMedianH->SetBinError(adcMedianH->FindBin(column),adcMeanUncert); 00664 adcSigmaH->SetBinContent(adcSigmaH->FindBin(column),adcSigma); 00665 adcSigmaH->SetBinError(adcSigmaH->FindBin(column),adcSigmaUncert); 00666 00667 meanAdcH->SetBinContent(ring,adcMean); 00668 meanAdcH->SetBinError(ring, adcMeanUncert); 00669 00670 gainChannelH->SetBinContent(ring*24+column,adcMean); 00671 gainChannelH->SetBinError(ring*24+column,adcMeanUncert); 00672 00673 00674 00675 } 00676 00677 // Record ring averages, ignoring bad channels 00678 if( occupancyOk && gainOk ){ 00679 if( occupancyUncert > 0 ){ 00680 meanOccRing[ring] += occupancy / (occupancyUncert*occupancyUncert); 00681 meanOccWeight[ring] += 1.0 / (occupancyUncert*occupancyUncert); 00682 } 00683 if( tdcOffsetUncert > 0 ){ 00684 meanTdcRing[ring] += tdcOffset / (tdcOffsetUncert*tdcOffsetUncert); 00685 meanTdcWeight[ring] += 1.0 / (tdcOffsetUncert*tdcOffsetUncert); 00686 } 00687 } 00688 00689 } 00690 00691 // Record mean time offset and mean occupancy by ring 00692 for(int ring = 1; ring <= nRings; ring++){ 00693 if( meanOccWeight[ring] > 0 ){ 00694 meanOccupancyH->SetBinContent(ring, 00695 meanOccRing[ring] / meanOccWeight[ring]); 00696 meanOccupancyH->SetBinError(ring, 00697 sqrt( 1.0 / meanOccWeight[ring] )); 00698 } 00699 if( meanTdcWeight[ring] > 0 ){ 00700 meanTdcOffsetH->SetBinContent(ring, 00701 meanTdcRing[ring] / meanTdcWeight[ring]); 00702 meanTdcOffsetH->SetBinError(ring, 00703 sqrt( 1.0 / meanTdcWeight[ring] )); 00704 } 00705 } // Loop over rings 00706 00708 //fill the data after cuts if the total events number less than max_event_number 00711 if(nReadouts<=max_event_number) 00712 for(int event_count=0;event_count<nReadouts;event_count++) 00713 { 00714 00715 //info()<<"test 2 event_count="<<event_count<<endreq; 00716 00717 TH1F* adcSumHWithCuts = m_statsSvc->getTH1F( this->getPath(detector) + "/adcSumWithCuts"); 00718 if(!adcSumHWithCuts) return StatusCode::FAILURE; 00719 TH1F* tdcMedianHWithCuts = m_statsSvc->getTH1F( this->getPath(detector)+"/tdcMedianWithCuts"); 00720 if(!tdcMedianHWithCuts) return StatusCode::FAILURE; 00721 00722 // Add the current readout to the calibration statistics 00723 std::vector<int> readoutTdcWithCuts; 00724 double adcSumWithCuts = 0; 00725 // Initialize statistics for each channel 00726 for(chanIter = channelList.begin(); chanIter != chanEnd; chanIter++) 00727 { 00728 DayaBay::FeeChannelId chanId = *chanIter; 00729 00730 // Analyze the channel data 00731 // Determine ring and column 00732 DayaBay::AdPmtSensor pmtId = m_cableSvc->adPmtSensor(chanId, svcMode); 00733 int ring = pmtId.ring(); 00734 int column = pmtId.column(); 00735 int channel_number=ring*100+column; 00736 00737 // Get Parameters and Histograms 00738 TObject* nHitsObjWithCuts = m_statsSvc->get(this->getPath(chanId) + "/nHitsWithCuts"); 00739 TParameter<int>* nHitsParWithCuts = dynamic_cast<TParameter<int>*>(nHitsObjWithCuts); 00740 if(!nHitsParWithCuts) return StatusCode::FAILURE; 00741 TH1F* tdcHWithCuts = m_statsSvc->getTH1F( this->getPath(chanId) + "/tdcWithCuts"); 00742 if(!tdcHWithCuts) return StatusCode::FAILURE; 00743 TH1F* adcRawHWithCuts = m_statsSvc->getTH1F( this->getPath(chanId) + "/adcRawWithCuts"); 00744 if(!adcRawHWithCuts) return StatusCode::FAILURE; 00745 TH1F* adcHWithCuts = m_statsSvc->getTH1F( this->getPath(chanId) + "/adcWithCuts"); 00746 if(!adcHWithCuts) return StatusCode::FAILURE; 00747 TH1F* adcByClockHWithCuts = m_statsSvc->getTH1F( this->getPath(chanId) 00748 + "/adcByClockWithCuts"); 00749 if(!adcByClockHWithCuts) return StatusCode::FAILURE; 00750 00751 // Loop over tdc values 00752 int chan_num=ring*24+column; 00753 double baseline=adc_baseline_buffer[chan_num%max_PMT_number][event_count]; 00754 for(int k=0;k<max_FEE_hit;k++) 00755 { 00756 int tdc=tdc_buffer[chan_num%max_PMT_number][event_count][k]; 00757 int adc=adc_buffer[chan_num%max_PMT_number][event_count][k]; 00758 int adcClock=adc_clock_buffer[chan_num%max_PMT_number][event_count][k]; 00759 bool hit_ok_flag=0; 00760 //if(tdc>0) 00761 //{ 00762 //info()<<"tdc="<<tdc<<",m_TDC_Cut_Min_Channel="<<m_TDC_Cut_Min_Channel[channel_number]<<",m_TDC_Cut_Max_Channel="<<m_TDC_Cut_Max_Channel[channel_number]<<endreq; 00763 //info()<<"adc="<<adc<<",m_ADC_Cut_Min_Channel="<<m_ADC_Cut_Min_Channel[channel_number]<<",m_ADC_Cut_Max_Channel="<<m_ADC_Cut_Max_Channel[channel_number]<<endreq; 00764 //} 00765 if(tdc>=m_TDC_Cut_Min_Channel[channel_number]&&tdc<=m_TDC_Cut_Max_Channel[channel_number]) 00766 if(adc>=m_ADC_Cut_Min_Channel[channel_number]&&adc<=m_ADC_Cut_Max_Channel[channel_number]) 00767 { 00768 //info()<<"test 5 chan_num="<<chan_num<<",event_count="<<event_count<<",k="<<k<<",adc="<<adc<<endreq; 00769 readoutTdcWithCuts.push_back(tdc); 00770 tdcHWithCuts->Fill(tdc); 00771 00772 adcRawHWithCuts->Fill(adc); 00773 adcHWithCuts->Fill(adc-baseline); 00774 adcByClockHWithCuts->Fill(adcClock,adc-baseline); 00775 adcSumWithCuts += adc-baseline; 00776 hit_ok_flag=1; 00777 00778 //info()<<"test 8 adcSumWithCuts="<<adcSumWithCuts<<endreq; 00779 } 00780 00781 00782 // Increment number of hits on this channel 00783 if(hit_ok_flag) 00784 nHitsParWithCuts->SetVal( nHitsParWithCuts->GetVal() + 1 ); 00785 } //over for for(int k=0;k<max_FEE_hit;k++) 00786 00787 //info()<<"test 9 adcSumWithCuts="<<adcSumWithCuts<<endreq; 00788 00789 00790 }// over for for(chanIter = channelList.begin(); chanIter != chanEnd; chanIter++) 00791 00792 if(readoutTdcWithCuts.size() > 0) 00793 { 00794 // Find the median TDC 00795 std::sort(readoutTdcWithCuts.begin(), readoutTdcWithCuts.end()); 00796 int medianIdx = readoutTdcWithCuts.size()/2; 00797 int medianTdc = readoutTdcWithCuts[medianIdx]; 00798 for(chanIter = channelList.begin(); chanIter != chanEnd; chanIter++) 00799 { 00800 DayaBay::FeeChannelId chanId = *chanIter; 00801 00802 // Analyze the channel data 00803 // Determine ring and column 00804 DayaBay::AdPmtSensor pmtId = m_cableSvc->adPmtSensor(chanId, svcMode); 00805 int ring = pmtId.ring(); 00806 int column = pmtId.column(); 00807 int channel_number=ring*100+column; 00808 00809 00810 TH1F* tdcByMedianHWithCuts = 00811 m_statsSvc->getTH1F( this->getPath(chanId) + "/tdcByMedianWithCuts"); 00812 if(!tdcByMedianHWithCuts) return StatusCode::FAILURE; 00813 00814 int chan_num=ring*24+column; 00815 for(int k=0;k<max_FEE_hit;k++) 00816 { 00817 int tdc=tdc_buffer[chan_num%max_PMT_number][event_count][k]; 00818 int adc=adc_buffer[chan_num%max_PMT_number][event_count][k]; 00819 00820 if(tdc>=m_TDC_Cut_Min_Channel[channel_number]&&tdc<=m_TDC_Cut_Max_Channel[channel_number]) 00821 if(adc>=m_ADC_Cut_Min_Channel[channel_number]&&adc<=m_ADC_Cut_Max_Channel[channel_number]) 00822 tdcByMedianHWithCuts->Fill(tdc-medianTdc); 00823 } 00824 } //over for for(chanIter = channelList.begin(); chanIter != chanEnd; chanIter++) 00825 tdcMedianHWithCuts->Fill(medianTdc); 00826 } //over for if(readoutTdcWithCuts.size() > 0) 00827 00828 //info()<<"test 3 event_count="<<event_count<<" adcSumWithCuts="<<adcSumWithCuts<<endreq; 00829 adcSumHWithCuts->Fill(adcSumWithCuts); 00830 00831 } //over if(nReadouts<=max_event_number) and for(int event_count=0;event_count<nReadouts;event_count++) 00832 00835 00836 00837 } // Loop over detectors 00838 00839 00840 return StatusCode::SUCCESS; 00841 }
| const InterfaceID & IPmtCalibParamTool::interfaceID | ( | ) | [static, inherited] |
Retrieve interface ID.
Reimplemented from IAlgTool.
Definition at line 8 of file IPmtCalibParamTool.cc.
00009 { 00010 return IID_IPmtCalibParamTool; 00011 }
std::string PmtCalibLeadingEdgeWithCuts::m_cableSvcName [private] |
Definition at line 128 of file PmtCalibLeadingEdgeWithCuts.h.
std::string PmtCalibLeadingEdgeWithCuts::m_calibSvcName [private] |
Definition at line 131 of file PmtCalibLeadingEdgeWithCuts.h.
std::string PmtCalibLeadingEdgeWithCuts::m_filepath [private] |
Definition at line 134 of file PmtCalibLeadingEdgeWithCuts.h.
ICableSvc* PmtCalibLeadingEdgeWithCuts::m_cableSvc [private] |
Definition at line 137 of file PmtCalibLeadingEdgeWithCuts.h.
Definition at line 140 of file PmtCalibLeadingEdgeWithCuts.h.
Definition at line 143 of file PmtCalibLeadingEdgeWithCuts.h.
std::vector<DayaBay::Detector> PmtCalibLeadingEdgeWithCuts::m_processedDetectors [private] |
Definition at line 146 of file PmtCalibLeadingEdgeWithCuts.h.
int PmtCalibLeadingEdgeWithCuts::adc_buffer[max_PMT_number][max_event_number][max_FEE_hit] [private] |
Definition at line 156 of file PmtCalibLeadingEdgeWithCuts.h.
int PmtCalibLeadingEdgeWithCuts::tdc_buffer[max_PMT_number][max_event_number][max_FEE_hit] [private] |
Definition at line 157 of file PmtCalibLeadingEdgeWithCuts.h.
int PmtCalibLeadingEdgeWithCuts::adc_clock_buffer[max_PMT_number][max_event_number][max_FEE_hit] [private] |
Definition at line 158 of file PmtCalibLeadingEdgeWithCuts.h.
double PmtCalibLeadingEdgeWithCuts::adc_baseline_buffer[max_PMT_number][max_event_number] [private] |
Definition at line 159 of file PmtCalibLeadingEdgeWithCuts.h.
double PmtCalibLeadingEdgeWithCuts::m_ADC_Cut_Min [private] |
Definition at line 166 of file PmtCalibLeadingEdgeWithCuts.h.
double PmtCalibLeadingEdgeWithCuts::m_ADC_Cut_Max [private] |
Definition at line 167 of file PmtCalibLeadingEdgeWithCuts.h.
std::map<int,double> PmtCalibLeadingEdgeWithCuts::m_ADC_Cut_Min_Channel [private] |
Definition at line 171 of file PmtCalibLeadingEdgeWithCuts.h.
std::map<int,double> PmtCalibLeadingEdgeWithCuts::m_ADC_Cut_Max_Channel [private] |
Definition at line 172 of file PmtCalibLeadingEdgeWithCuts.h.
double PmtCalibLeadingEdgeWithCuts::m_TDC_Cut_Min [private] |
Definition at line 179 of file PmtCalibLeadingEdgeWithCuts.h.
double PmtCalibLeadingEdgeWithCuts::m_TDC_Cut_Max [private] |
Definition at line 180 of file PmtCalibLeadingEdgeWithCuts.h.
std::map<int,double> PmtCalibLeadingEdgeWithCuts::m_TDC_Cut_Min_Channel [private] |
Definition at line 183 of file PmtCalibLeadingEdgeWithCuts.h.
std::map<int,double> PmtCalibLeadingEdgeWithCuts::m_TDC_Cut_Max_Channel [private] |
Definition at line 184 of file PmtCalibLeadingEdgeWithCuts.h.
1.4.7