/search.css" rel="stylesheet" type="text/css"/> /search.js">
| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |

In This Package:

RandomGainAlg.cc
Go to the documentation of this file.
00001 #include "RandomGainAlg.h"
00002 #include "NIMmodel.h"
00003 #include "Event/ReadoutHeader.h"
00004 #include "Event/Readout.h"
00005 #include "Event/ReadoutPmtCrate.h"
00006 #include "Event/ReadoutRpcCrate.h"
00007 #include "TH1F.h"
00008 #include "TF1.h"
00009 #include "TH2.h"
00010 #include "TParameter.h"
00011 #include <math.h>
00012 #include "TGraph.h"
00013 #include "TMath.h"
00014 
00015 #include "CLHEP/Units/SystemOfUnits.h"
00016 
00017 using namespace DayaBay;
00018 using namespace std;
00019 
00020 
00021 RandomGainAlg::RandomGainAlg(const std::string& name, ISvcLocator* pSvcLocator)
00022     : GaudiAlgorithm(name,pSvcLocator)
00023 {
00024         declareProperty("ReadoutLocation",m_readoutLocation=DayaBay::ReadoutHeaderLocation::Default,"Location in the TES where the input ReadoutHeader is to be found.");
00025         declareProperty("CableSvcName",m_cableSvcName="CableSvc", "Name of service to map between detector, hardware, and electronic IDs");
00026         declareProperty("FilePath",m_filepath="/file0/pmtCalibRandomGain", "File path of registered histograms.");
00027 }
00028 
00029 RandomGainAlg::~RandomGainAlg(){}
00030 
00031 StatusCode RandomGainAlg::initialize(){
00032         
00033         info() << "initialize()" << endreq;
00034         m_cableSvc = svc<ICableSvc>(m_cableSvcName, true);
00035         if(!m_cableSvc){
00036         error() << "Failed to access cable svc: " << m_cableSvcName << endreq;
00037         return StatusCode::FAILURE;
00038         }
00039   
00040         m_statsSvc = svc<IStatisticsSvc>("StatisticsSvc",true);
00041         if(!m_statsSvc) {
00042         error() << " No StatisticsSvc available." << endreq;
00043         return StatusCode::FAILURE;
00044         }
00045 
00046         std::vector<std::string> detDirs = m_statsSvc->getSubFolders(m_filepath);
00047         std::vector<std::string>::iterator dirIter, dirEnd = detDirs.end();
00048         for(dirIter=detDirs.begin(); dirIter != dirEnd; dirIter++){
00049         info() << "Processing input directory: " << m_filepath << "/"  << *dirIter << endreq;
00050         DayaBay::Detector detector( DayaBay::Detector::siteDetPackedFromString(*dirIter) );
00051         if( detector.site() == Site::kUnknown || detector.detectorId() == DetectorId::kUnknown ){
00052                 warning() << "Input directory: " << m_filepath << "/" << *dirIter << " not processed." << endreq;
00053         }
00054                 m_processedDetectors.push_back(detector);
00055         }
00056         m_first = true;
00057         m_lastTriTimeAD1=0;
00058         m_lastTriTimeAD2=0;
00059         m_lastTriTimeWPI=0;
00060         m_lastTriTimeWPO=0;
00061         total_event=0;
00062         used_event =0;
00063         return StatusCode::SUCCESS;
00064 }
00065 
00066 StatusCode RandomGainAlg::execute(){
00067         debug() << "calling execute()" << endreq; 
00068         DayaBay::ReadoutHeader* readoutHeader = get<DayaBay::ReadoutHeader>(m_readoutLocation);
00069         TimeStamp trigTime = readoutHeader->timeStamp();
00070 
00071         if(m_first) {
00072                 m_beginTime = trigTime;
00073         m_first = false;
00074         }
00075         m_endTime = trigTime;
00076         
00077         DayaBay::Detector det(readoutHeader->context().GetSite(),readoutHeader->context().GetDetId());
00078         debug() << "Got readout from " << det.detName() <<" id= "<< det.siteDetPackedData() <<endreq;
00079         
00080         if(det.detectorId()==7) return StatusCode::SUCCESS;
00081         
00082         if(det.isAD() || det.isWaterShield()){
00083                 const DayaBay::DaqCrate* daqCrate = readoutHeader->daqCrate();
00084                 if(!daqCrate){
00085                 error() << "Failed to get DAQ crate from header" << endreq;
00086                 return StatusCode::FAILURE;
00087         }
00088         const DayaBay::DaqPmtCrate* pmtCrate= daqCrate->asPmtCrate();
00089                 if(!pmtCrate){
00090                 debug() << "process(): Do not process detector: " << det.detName() << endreq;
00091                 return StatusCode::SUCCESS;
00092                 }
00093                 int trigType = pmtCrate->triggerType();
00094                 if(det.detectorId()==1&&trigType == 0x10000020) total_event++;
00095                 if(det.detectorId()==1&&(trigTime.GetNanoSec()-m_lastTriTimeAD1)<20000) { m_lastTriTimeAD1 = trigTime.GetNanoSec(); return StatusCode::SUCCESS;}
00096                 if(det.detectorId()==2&&(trigTime.GetNanoSec()-m_lastTriTimeAD2)<20000) { m_lastTriTimeAD2 = trigTime.GetNanoSec(); return StatusCode::SUCCESS;}
00097                 if(det.detectorId()==5&&(trigTime.GetNanoSec()-m_lastTriTimeWPI)<20000) { m_lastTriTimeWPI = trigTime.GetNanoSec(); return StatusCode::SUCCESS;}
00098                 if(det.detectorId()==6&&(trigTime.GetNanoSec()-m_lastTriTimeWPO)<20000) { m_lastTriTimeWPO = trigTime.GetNanoSec(); return StatusCode::SUCCESS;}
00099         
00100                 if(det.detectorId()==1) { m_lastTriTimeAD1 = trigTime.GetNanoSec();}
00101                 if(det.detectorId()==2) { m_lastTriTimeAD2 = trigTime.GetNanoSec();}
00102                 if(det.detectorId()==5) { m_lastTriTimeWPI = trigTime.GetNanoSec(); }
00103                 if(det.detectorId()==6) { m_lastTriTimeWPO = trigTime.GetNanoSec();}
00104         
00106                 
00107                 if(trigType != 0x10000020) return StatusCode::SUCCESS;
00108                 
00109                 Context context = readoutHeader->context();
00110                 if( !hasStats(det) ){
00111                         StatusCode sc = prepareStats(context);
00112                 if( sc != StatusCode::SUCCESS ) return sc;
00113                 }
00114                 ServiceMode svcMode(context, 0);
00115                 const DayaBay::DaqPmtCrate::PmtChannelPtrList& channels= pmtCrate->channelReadouts();
00116         DayaBay::DaqPmtCrate::PmtChannelPtrList::const_iterator channelIter, channelEnd = channels.end();
00117         if(channels.size()>15) return StatusCode::SUCCESS;
00118         used_event++;
00119                 for(channelIter = channels.begin(); channelIter!=channelEnd; channelIter++) { 
00120                 const DayaBay::DaqPmtChannel& channel = *(*channelIter);
00121                         const DayaBay::FeeChannelId& chanId = channel.channelId();
00122                         const DayaBay::DetectorSensor& sensDetId 
00123                 = m_cableSvc->sensor(chanId, svcMode);
00124                 if (sensDetId.bogus()) {
00125                  warning() << "Got bogus PMT: \"" << sensDetId
00126                           << "\" using channelId: \"" << chanId
00127                           << "\" and context \"" << context.AsString() << "\""
00128                           << endreq;
00129                     continue;
00130                 }
00131 
00132                         //if((channel.tdc()).size()>1) continue;
00133                         if(channel.hitCount()>1) continue;
00134                         TH1F* tdcRawH = m_statsSvc->getTH1F( this->getPath(chanId) + "/tdcRaw");
00135                 if(!tdcRawH) return StatusCode::FAILURE;
00136                 TH1F* adcRawH = m_statsSvc->getTH1F( this->getPath(chanId) + "/adcRaw");
00137                 if(!adcRawH) return StatusCode::FAILURE;
00138                 TH1F* preAdcH = m_statsSvc->getTH1F( this->getPath(chanId) + "/preAdc");
00139                 if(!preAdcH) return StatusCode::FAILURE;
00140                 TH1F* adcH = m_statsSvc->getTH1F( this->getPath(chanId) + "/adc");    
00141                 if(!adcH) return StatusCode::FAILURE;
00142                         
00143                         int tdc = channel.tdc(0);
00144                         int adc = channel.adc(0);
00145                         float preAdc = channel.preAdcAvg(0);
00146                         int cycle = channel.peakCycle(0);
00147                         
00148                         if(channel.isHighGainAdc(0)&&cycle>3&&cycle<8 ) {
00149                                 tdcRawH->Fill(tdc);
00150                                 adcRawH->Fill(adc);
00151                                 preAdcH->Fill(preAdc);
00152                                 adcH->Fill(adc-preAdc);
00153                         }
00154                 }
00155         }
00156         else
00157                 return StatusCode::SUCCESS;
00158                 
00159         return StatusCode::SUCCESS;             
00160 }
00161 
00162 StatusCode RandomGainAlg::finalize(){
00163         info()<<"Begin time: "<< m_beginTime.GetSec()<< endreq; 
00164         info()<<"End Time "<< m_endTime.GetSec() <<endreq;
00165         info()<<"Total Random Triggers in AD1" << total_event<<endreq;
00166         info()<<"used Random Triggers in AD1" << used_event<<endreq;
00167         
00168         std::vector<DayaBay::Detector>::iterator detIter, detEnd = m_processedDetectors.end();
00169         for(detIter = m_processedDetectors.begin(); detIter != detEnd; detIter++){
00170         DayaBay::Detector detector = *detIter;
00171         if(detector.isAD()) {
00172                 TH1F* expAdcMeanH = m_statsSvc->getTH1F( this->getPath(detector)+"/expAdcMean");
00173                 if(!expAdcMeanH) return StatusCode::FAILURE;
00174                 TH1F* expAdcSigmaH = m_statsSvc->getTH1F( this->getPath(detector)+"/expAdcSigma");
00175                 if(!expAdcSigmaH) return StatusCode::FAILURE;
00176                 TH2F* gainH2D = m_statsSvc->getTH2F( this->getPath(detector) + "/ADgain2D");
00177                 if(!gainH2D) return StatusCode::FAILURE;
00178                 TH1F* gainH1D = m_statsSvc->getTH1F( this->getPath(detector)  + "/gain1D");
00179                 if(!gainH1D) return StatusCode::FAILURE;
00180                         TH2F* Chi2NDF = m_statsSvc->getTH2F( this->getPath(detector) + "/ADChi2NDF");
00181                 if(!Chi2NDF) return StatusCode::FAILURE;
00182                 TH1F* Chi2NDF_1D = m_statsSvc->getTH1F( this->getPath(detector) + "/Chi2_1D");
00183                 if(!Chi2NDF_1D) return StatusCode::FAILURE;
00184                         
00185                         TObject* simFlagObj = m_statsSvc->get(this->getPath(detector) +"/simFlag");
00186                 TParameter<int>* simFlagPar = dynamic_cast<TParameter<int>*>(simFlagObj);
00187                 if(!simFlagPar) return StatusCode::FAILURE;
00188 
00189                 TObject* timeSecObj = m_statsSvc->get(this->getPath(detector) +"/timeSec");
00190                 TParameter<int>* timeSecPar = dynamic_cast<TParameter<int>*>(timeSecObj);
00191                 if(!timeSecPar) return StatusCode::FAILURE;
00192 
00193                 TObject* timeNanoSecObj = 
00194                  m_statsSvc->get(this->getPath(detector) +"/timeNanoSec");
00195                 TParameter<int>* timeNanoSecPar = 
00196                 dynamic_cast<TParameter<int>*>(timeNanoSecObj);
00197                 if(!timeNanoSecPar) return StatusCode::FAILURE;
00198 
00199                 Site::Site_t site = detector.site();
00200                         DetectorId::DetectorId_t detId = detector.detectorId();
00201                 SimFlag::SimFlag_t simFlag = (SimFlag::SimFlag_t)(simFlagPar->GetVal());
00202                 time_t timeSec = (time_t)(timeSecPar->GetVal());
00203                 int timeNanoSec = timeNanoSecPar->GetVal();
00204 
00205                 Context context(site,simFlag,TimeStamp(timeSec,timeNanoSec),detId);
00206                 ServiceMode svcMode(context, 0);
00207 
00208                         const std::vector<DayaBay::FeeChannelId>& channelList = m_cableSvc->feeChannelIds( svcMode );
00209                 std::vector<DayaBay::FeeChannelId>::const_iterator chanIter, chanEnd = channelList.end();
00210 
00211                     for(chanIter = channelList.begin(); chanIter != chanEnd; chanIter++){
00212                         DayaBay::FeeChannelId chanId = *chanIter;
00213                         TH1F* adcH = m_statsSvc->getTH1F( this->getPath(chanId) + "/adc");
00214                         if(!adcH) return StatusCode::FAILURE;
00215                         DayaBay::AdPmtSensor pmtId = m_cableSvc->adPmtSensor(chanId, svcMode);
00216                         int ring = pmtId.ring();
00217                         int column = pmtId.column();
00218                         double adcMean = 0;
00219                         double adcMeanUncert = 0;
00220                         double adcSigma = 0;
00221                         double adcSigmaUncert = 0;
00222                                 TF1 *nimmodelfunc = new TF1("nimmodelfunc",NIMmodel,5,100,8);
00223         nimmodelfunc->SetParNames(  "N","Q_{0}","#sigma_{0}","Q_{1}","#sigma_{1}","w", "#alpha","#mu");
00224         double AreaGuess, Q1MeanGuess, Q1SigmaGuess;
00225         float q0=0;//for a pedestal substracted distribution
00226         float sigma0=0;//for a pedestal substracted distribution
00227         AreaGuess = adcH->GetEntries();
00228         Q1MeanGuess = adcH->GetMean() - q0;
00229         Q1SigmaGuess = adcH->GetRMS();
00230         nimmodelfunc->SetParameters(AreaGuess, 90,1.9,Q1MeanGuess, Q1SigmaGuess, 0.01, 0.03,    0.1);
00231         //fix four parameters
00232         nimmodelfunc->FixParameter(1,q0);
00233         nimmodelfunc->FixParameter(2,sigma0);
00234         nimmodelfunc->FixParameter(5,0.);
00235         nimmodelfunc->FixParameter(6,0.);
00236        
00237         adcH->Fit(nimmodelfunc,"R");
00238         if(nimmodelfunc!=NULL){
00239           
00240           adcMean=nimmodelfunc->GetParameter(3);
00241           adcMeanUncert =nimmodelfunc->GetParError(3);
00242           adcSigma = nimmodelfunc->GetParameter(4);
00243           adcSigmaUncert = nimmodelfunc->GetParError(4);
00244           
00245           if(fabs(adcH->GetMean()-nimmodelfunc->GetParameter(3)) > adcH->GetRMS()){
00246             warning() << "Issues with fit on ring, column: " << ring << " / " 
00247                       << column << " mean/peak: " << adcMean << " / " << adcH->GetMean() 
00248                       << endreq;
00249             //gainOk = false;
00250           }//print warning if issues
00251           
00252           
00253         } else {
00254 
00255           adcMean=0;
00256           adcMeanUncert=0;
00257           adcSigma = 0;
00258           adcSigmaUncert = 0;
00259 
00260           warning() << "PMT full model fit error on ring/column: " << ring << "," << column << endreq;
00261           //gainOk = false;
00262         }
00263          double Chi2 = nimmodelfunc->GetChisquare();
00264         double ndf = nimmodelfunc->GetNDF();
00265         double fitquality = Chi2/ndf;
00266                                 
00267                                 
00268         
00269                                 gainH2D->SetBinContent(gainH2D->GetXaxis()->FindBin(column),
00270                                         gainH2D->GetYaxis()->FindBin(ring),
00271                                         adcMean);
00272                                 gainH2D->SetBinError(gainH2D->GetXaxis()->FindBin(column),
00273                                         gainH2D->GetYaxis()->FindBin(ring),
00274                                         adcMeanUncert);
00275                                 gainH1D->Fill(adcMean);
00276                                 Chi2NDF->SetBinContent(Chi2NDF->GetXaxis()->FindBin(column),
00277                                         Chi2NDF->GetYaxis()->FindBin(ring),
00278                                         fitquality);
00279                                 if(chanId.board()<=18) {
00280                                         expAdcMeanH->SetBinContent(chanId.board()*16+chanId.connector(),adcMean);
00281                                         expAdcMeanH->SetBinError(chanId.board()*16+chanId.connector(),adcMeanUncert);
00282                                         expAdcSigmaH->SetBinContent(chanId.board()*16+chanId.connector(),adcSigma);
00283                                         expAdcSigmaH->SetBinError(chanId.board()*16+chanId.connector(),adcSigmaUncert);
00284                                         Chi2NDF_1D->SetBinContent(chanId.board()*16+chanId.connector(), fitquality);
00285                                 }
00286                         }
00287                 }
00288                 else if(detector.isWaterShield()){  
00289                         TH1F* expAdcMeanH = m_statsSvc->getTH1F( this->getPath(detector)+"/expAdcMean");
00290                 if(!expAdcMeanH) return StatusCode::FAILURE;
00291                 TH1F* expAdcSigmaH = m_statsSvc->getTH1F( this->getPath(detector)+"/expAdcSigma");
00292                 if(!expAdcSigmaH) return StatusCode::FAILURE;
00293                 TH1F* gainH1D = m_statsSvc->getTH1F( this->getPath(detector)  + "/gain1D");
00294                 if(!gainH1D) return StatusCode::FAILURE;
00295                 TH1F* Chi2NDF_1D = m_statsSvc->getTH1F( this->getPath(detector) + "/Chi2_1D");
00296                 if(!Chi2NDF_1D) return StatusCode::FAILURE;
00297                 
00298                 TObject* simFlagObj = m_statsSvc->get(this->getPath(detector) +"/simFlag");
00299                 TParameter<int>* simFlagPar = dynamic_cast<TParameter<int>*>(simFlagObj);
00300                 if(!simFlagPar) return StatusCode::FAILURE;
00301 
00302                 TObject* timeSecObj = m_statsSvc->get(this->getPath(detector) +"/timeSec");
00303                 TParameter<int>* timeSecPar = dynamic_cast<TParameter<int>*>(timeSecObj);
00304                 if(!timeSecPar) return StatusCode::FAILURE;
00305 
00306                 TObject* timeNanoSecObj = 
00307                  m_statsSvc->get(this->getPath(detector) +"/timeNanoSec");
00308                 TParameter<int>* timeNanoSecPar = 
00309                 dynamic_cast<TParameter<int>*>(timeNanoSecObj);
00310                 if(!timeNanoSecPar) return StatusCode::FAILURE;
00311 
00312                 Site::Site_t site = detector.site();
00313                         DetectorId::DetectorId_t detId = detector.detectorId();
00314                 SimFlag::SimFlag_t simFlag = (SimFlag::SimFlag_t)(simFlagPar->GetVal());
00315                 time_t timeSec = (time_t)(timeSecPar->GetVal());
00316                 int timeNanoSec = timeNanoSecPar->GetVal();
00317 
00318                 Context context(site,simFlag,TimeStamp(timeSec,timeNanoSec),detId);
00319                 ServiceMode svcMode(context, 0);
00320                 
00321                 const std::vector<DayaBay::FeeChannelId>& channelList = m_cableSvc->feeChannelIds( svcMode );
00322                 std::vector<DayaBay::FeeChannelId>::const_iterator chanIter, chanEnd = channelList.end();
00323 
00324                     for(chanIter = channelList.begin(); chanIter != chanEnd; chanIter++){
00325                         DayaBay::FeeChannelId chanId = *chanIter;
00326                         TH1F* adcH = m_statsSvc->getTH1F( this->getPath(chanId) + "/adc");
00327                         if(!adcH) return StatusCode::FAILURE;
00328                         //DayaBay::AdPmtSensor pmtId = m_cableSvc->adPmtSensor(chanId, svcMode);
00329                         //int ring = pmtId.ring();
00330                         //int column = pmtId.column();
00331                         double adcMean = 0;
00332                         double adcMeanUncert = 0;
00333                         double adcSigma = 0;
00334                         double adcSigmaUncert = 0;
00335 
00336                                 TF1 *nimmodelfunc = new TF1("nimmodelfunc",NIMmodel,5,100,8);
00337                                 nimmodelfunc->SetParNames(  "N","Q_{0}","#sigma_{0}","Q_{1}","#sigma_{1}","w", "#alpha","#mu");
00338                                 double AreaGuess, Q1MeanGuess, Q1SigmaGuess;
00339                                 float q0=0;//for a pedestal substracted distribution
00340                                 float sigma0=0;//for a pedestal substracted distribution
00341                                 AreaGuess = adcH->GetEntries();
00342                                 Q1MeanGuess = adcH->GetMean() - q0;
00343                                 Q1SigmaGuess = adcH->GetRMS();
00344                                 nimmodelfunc->SetParameters(AreaGuess, 90,1.9,Q1MeanGuess, Q1SigmaGuess, 0.01, 0.03,    0.1);
00345                                 //fix four parameters
00346                                 nimmodelfunc->FixParameter(1,q0);
00347                                 nimmodelfunc->FixParameter(2,sigma0);
00348                                 nimmodelfunc->FixParameter(5,0.);
00349                                 nimmodelfunc->FixParameter(6,0.);
00350        
00351                                 adcH->Fit(nimmodelfunc,"R");
00352                                 if(nimmodelfunc!=NULL){
00353                                  adcMean=nimmodelfunc->GetParameter(3);
00354                                  adcMeanUncert =nimmodelfunc->GetParError(3);
00355                                  adcSigma = nimmodelfunc->GetParameter(4);
00356                                  adcSigmaUncert = nimmodelfunc->GetParError(4);
00357                                  if(fabs(adcH->GetMean()-nimmodelfunc->GetParameter(3)) > adcH->GetRMS()){
00358                                 //warning() << "Issues with fit on ring, column: " << ring << " / " 
00359                                 //<< column << " mean/peak: " << adcMean << " / " << adcH->GetMean() 
00360                                 //<< endreq;
00361                                 }//print warning if issues
00362                           } else {
00363                                         adcMean=0;
00364                                         adcMeanUncert=0;
00365                                         adcSigma = 0;
00366                                         adcSigmaUncert = 0;
00367 
00368                                         //warning() << "PMT full model fit error on ring/column: " << ring << "," << column << endreq;
00369                                 }
00370                                 double Chi2 = nimmodelfunc->GetChisquare();
00371                                 double ndf = nimmodelfunc->GetNDF();
00372                                 double fitquality = Chi2/ndf;
00373 
00374                                 
00375                                 gainH1D->Fill(adcMean);
00376                                 if(chanId.board()<=18) {
00377                                         expAdcMeanH->SetBinContent(chanId.board()*16+chanId.connector(),adcMean);
00378                                         expAdcMeanH->SetBinError(chanId.board()*16+chanId.connector(),adcMeanUncert);
00379                                         expAdcSigmaH->SetBinContent(chanId.board()*16+chanId.connector(),adcSigma);
00380                                         expAdcSigmaH->SetBinError(chanId.board()*16+chanId.connector(),adcSigmaUncert);
00381                                         Chi2NDF_1D->SetBinContent(chanId.board()*16+chanId.connector(), fitquality);
00382                                 }
00383                         }
00384                 }
00385                 else 
00386                         return StatusCode::SUCCESS;
00387         }//end of detectors
00388         
00389         return StatusCode::SUCCESS;
00390 }
00391                 
00392 bool RandomGainAlg::hasStats(const DayaBay::Detector& detector){
00393   // Check if statistics have been initialized for this detector
00394   return (std::find(m_processedDetectors.begin(), 
00395                     m_processedDetectors.end(), 
00396                     detector) 
00397           != m_processedDetectors.end());
00398 }
00399 
00400 std::string RandomGainAlg::getPath( const DayaBay::Detector& detector ){
00401   return m_filepath + "/" + detector.detName();
00402 }
00403 
00404 std::string RandomGainAlg::getPath(const DayaBay::FeeChannelId& chanId){
00405   std::ostringstream path;
00406   path << m_filepath << "/" << chanId.detName() 
00407        << "/chan_" << chanId.board() << "_" << chanId.connector(); 
00408   return path.str();
00409 }
00410 
00411 StatusCode RandomGainAlg::prepareStats(const Context& context){
00412 
00413         DayaBay::Detector detector(context.GetSite(), context.GetDetId());
00414         // SimFlag
00415         {
00416         std::string name = "simFlag";
00417         std::ostringstream path;
00418         path << this->getPath(detector) << "/" << name;
00419         TParameter<int>* par = new TParameter<int>(name.c_str(), 
00420                                                context.GetSimFlag());
00421         if( m_statsSvc->put(path.str(),par).isFailure() ) {
00422                 error() << "prepareStats(): Could not register " << name 
00423                         << " at " << path << endreq;
00424                 delete par;
00425                 par = 0;
00426                 return StatusCode::FAILURE;
00427         }
00428         }
00429         
00430         // Timestamp: seconds
00431         {
00432         std::string name = "timeSec";
00433         std::ostringstream path;
00434         path << this->getPath(detector) << "/" << name;
00435         TParameter<int>* par = new TParameter<int>(name.c_str(), 
00436                                                context.GetTimeStamp().GetSec());
00437         if( m_statsSvc->put(path.str(),par).isFailure() ) {
00438                 error() << "prepareStats(): Could not register " << name 
00439                         << " at " << path << endreq;
00440                 delete par;
00441                 par = 0;
00442                 return StatusCode::FAILURE;
00443         }
00444         }
00445 
00446   // Timestamp: nanoseconds
00447         {
00448         std::string name = "timeNanoSec";
00449         std::ostringstream path;
00450         path << this->getPath(detector) << "/" << name;
00451         TParameter<int>* par = new TParameter<int>(name.c_str(), context.GetTimeStamp().GetNanoSec());
00452         if( m_statsSvc->put(path.str(),par).isFailure() ) {
00453                 error() << "prepareStats(): Could not register " << name 
00454                         << " at " << path << endreq;
00455                 delete par;
00456                 par = 0;
00457                 return StatusCode::FAILURE;
00458         }
00459         }
00460         
00461         // Adc - Ave preAdc, mean
00462   {
00463     std::ostringstream title, path;
00464     std::string name = "expAdcMean";
00465     title << "ADC - ave preAdc, mean in " << detector.detName();
00466     path << this->getPath(detector) << "/" << name;
00467     TH1F* expAdcMean = new TH1F(name.c_str(),title.str().c_str(),400,0,400); // to hold water PMT fee board 
00468     expAdcMean->GetXaxis()->SetTitle("board*16+connector");
00469     expAdcMean->GetYaxis()->SetTitle("ExpAdc mean");
00470     if( m_statsSvc->put(path.str(),expAdcMean).isFailure() ) {
00471       error() << "prepareStats(): Could not register " << path << endreq;
00472       delete expAdcMean;
00473       return StatusCode::FAILURE;
00474     }
00475   }
00476   // Adc - Ave preAdc, sigma
00477   {
00478     std::ostringstream title, path;
00479     std::string name = "expAdcSigma";
00480     title << "ADC - ave preAdc, sigma in " << detector.detName();
00481     path << this->getPath(detector) << "/" << name;
00482     TH1F* expAdcSigma = new TH1F(name.c_str(),title.str().c_str(),400,0,400);
00483     expAdcSigma->GetXaxis()->SetTitle("board*16+connector");
00484     expAdcSigma->GetYaxis()->SetTitle("ExpAdc sigma");
00485     if( m_statsSvc->put(path.str(),expAdcSigma).isFailure() ) {
00486       error() << "prepareStats(): Could not register " << path << endreq;
00487       delete expAdcSigma;
00488       return StatusCode::FAILURE;
00489     }
00490   }
00491     // fitting chi2/ndf
00492   {
00493     std::ostringstream title, path;
00494     std::string name = "Chi2_1D";
00495     title << "Chi2/NDF of fitting " << detector.detName();
00496     path << this->getPath(detector) << "/" << name;
00497     TH1F* Chi2_1D = new TH1F(name.c_str(),title.str().c_str(),400,0,400);
00498     Chi2_1D->GetXaxis()->SetTitle("board*16+connector");
00499     Chi2_1D->GetYaxis()->SetTitle("Chi2/NDF");
00500     if( m_statsSvc->put(path.str(),Chi2_1D).isFailure() ) {
00501       error() << "prepareStats(): Could not register " << path << endreq;
00502       delete Chi2_1D;
00503       return StatusCode::FAILURE;
00504     }
00505   }
00506 
00507         // AD Gain for each channel (jpochoa)
00508   {
00509       std::ostringstream title, path;
00510       std::string name = "ADgain2D";
00511       title << "AD Gain per channel " << detector.detName();
00512       path << this->getPath(detector) << "/" << name;
00513       TH2F* ADgain2D = new TH2F(name.c_str(),title.str().c_str(),49,0.25,24.75,19,-0.75,8.75);
00514       ADgain2D->GetXaxis()->SetTitle("Column");
00515       ADgain2D->GetYaxis()->SetTitle("Ring");
00516       if( m_statsSvc->put(path.str(),ADgain2D).isFailure() ) {
00517         error() << "prepareStats(): Could not register " << path << endreq;
00518         delete ADgain2D;
00519         return StatusCode::FAILURE;
00520       }
00521     }
00522 
00523         // Gain for all channels (jpochoa)
00524   {
00525       std::ostringstream title, path;
00526       std::string name = "gain1D";
00527       title << "Gain for all channels " << detector.detName();
00528       path << this->getPath(detector) << "/" << name;
00529       TH1F* gain1D = new TH1F(name.c_str(),title.str().c_str(),400,-200,200);
00530       gain1D->GetXaxis()->SetTitle("ADC/PE");
00531       gain1D->GetYaxis()->SetTitle("Number of entries");
00532       if( m_statsSvc->put(path.str(),gain1D).isFailure() ) {
00533         error() << "prepareStats(): Could not register " << path << endreq;
00534         delete gain1D;
00535         return StatusCode::FAILURE;
00536       }
00537     }
00538     
00539     // AD Chi2/ndf for each channel (yuzy)
00540   {
00541       std::ostringstream title, path;
00542       std::string name = "ADChi2NDF";
00543       title << "AD Chi2/NDF per channel " << detector.detName();
00544       path << this->getPath(detector) << "/" << name;
00545       TH2F* ADChi2NDF = new TH2F(name.c_str(),title.str().c_str(),49,0.25,24.75,19,-0.75,8.75);
00546       ADChi2NDF->GetXaxis()->SetTitle("Column");
00547       ADChi2NDF->GetYaxis()->SetTitle("Ring");
00548       if( m_statsSvc->put(path.str(),ADChi2NDF).isFailure() ) {
00549         error() << "prepareStats(): Could not register " << path << endreq;
00550         delete ADChi2NDF;
00551         return StatusCode::FAILURE;
00552       }
00553     }
00554 
00555         ServiceMode svcMode(context, 0);
00556 
00557   // Get list of all FEE channels
00558   const std::vector<DayaBay::FeeChannelId>& channelList 
00559     = m_cableSvc->feeChannelIds( svcMode );
00560   std::vector<DayaBay::FeeChannelId>::const_iterator chanIter, 
00561     chanEnd = channelList.end();
00562   // Initialize statistics for each channel
00563   for(chanIter = channelList.begin(); chanIter != chanEnd; chanIter++){
00564     DayaBay::FeeChannelId chanId = *chanIter;
00565         // Raw TDC spectrum
00566     {
00567       std::string name = "tdcRaw";
00568       std::ostringstream title, path;
00569       title << "Raw TDC values for " << detector.detName() 
00570             << " channel " << chanId.board() << "_"
00571             << chanId.connector();
00572       path << this->getPath(chanId) << "/" << name;
00573       TH1F* tdcRaw = new TH1F(name.c_str(),title.str().c_str(),
00574                               2000,0,2000);
00575       tdcRaw->GetXaxis()->SetTitle("TDC value");
00576       tdcRaw->GetYaxis()->SetTitle("Entries");
00577       if( m_statsSvc->put(path.str(),tdcRaw).isFailure() ) {
00578         error() << "prepareStats(): Could not register " << path << endreq;
00579         delete tdcRaw;
00580         return StatusCode::FAILURE;
00581       }
00582     }
00583         // Raw ADC spectrum
00584     {
00585       std::ostringstream title, path;
00586       std::string name = "adcRaw";
00587       title << "ADC values for " << detector.detName() 
00588             << " channel " << chanId.board() << "_"
00589             << chanId.connector();
00590       path << this->getPath(chanId) << "/" << name;
00591       TH1F* adcRaw = new TH1F(name.c_str(),title.str().c_str(),
00592                               4096,0,4096);
00593       adcRaw->GetXaxis()->SetTitle("ADC value");
00594       adcRaw->GetYaxis()->SetTitle("Entries");
00595       if( m_statsSvc->put(path.str(),adcRaw).isFailure() ) {
00596         error() << "prepareStats(): Could not register " << path << endreq;
00597         delete adcRaw;
00598         return StatusCode::FAILURE;
00599       }
00600     }
00601 
00602         // ADC spectrum, average pedestal substracted
00603     {
00604       std::ostringstream title, path;
00605       std::string name = "adc";
00606       title << "ADC - Average preAdc " << detector.detName() 
00607             << " channel " << chanId.board() << "_"
00608             << chanId.connector();
00609       path << this->getPath(chanId) << "/" << name;
00610       TH1F* adc = new TH1F(name.c_str(),title.str().c_str(),
00611                            120,-20,100);
00612       adc->GetXaxis()->SetTitle("ADC value");
00613       adc->GetYaxis()->SetTitle("Entries");
00614       if( m_statsSvc->put(path.str(),adc).isFailure() ) {
00615         error() << "prepareStats(): Could not register " << path << endreq;
00616         delete adc;
00617         return StatusCode::FAILURE;
00618       }
00619     }
00620 
00621     // PreADC spectrum, baseline subtracted
00622     {
00623       std::ostringstream title, path;
00624       std::string name = "preAdc";
00625       title << "PreADC values for " << detector.detName()
00626             << " channel " << chanId.board() << "_"
00627             << chanId.connector();
00628       path << this->getPath(chanId) << "/" << name;
00629       TH1F* preAdc = new TH1F(name.c_str(),title.str().c_str(),
00630                            1200,-200,1000);
00631       preAdc->GetXaxis()->SetTitle("PreADC value");
00632       preAdc->GetYaxis()->SetTitle("Entries");
00633       if( m_statsSvc->put(path.str(),preAdc).isFailure() ) {
00634         error() << "prepareStats(): Could not register " << path << endreq;
00635         delete preAdc;
00636         return StatusCode::FAILURE;
00637       }
00638     }
00639         }
00640         
00641         m_processedDetectors.push_back(detector);
00642 
00643   return StatusCode::SUCCESS;
00644 }
00645
| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |
Generated on Fri May 16 2014 10:03:04 for RandomGain by doxygen 1.7.4