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

In This Package:

Public Member Functions | Protected Member Functions | Protected Attributes | Private Attributes
CalibStatsAlg Class Reference

#include <CalibStatsAlg.h>

Inheritance diagram for CalibStatsAlg:
Inheritance graph
[legend]
Collaboration diagram for CalibStatsAlg:
Collaboration graph
[legend]

List of all members.

Public Member Functions

 CalibStatsAlg (const std::string &name, ISvcLocator *pSvcLocator)
virtual ~CalibStatsAlg ()
virtual StatusCode initialize ()
virtual StatusCode execute ()
virtual StatusCode finalize ()
DayaBay::UserDataHeaderGetCurrentHeaderObject () const
virtual StatusCode sysInitialize ()
virtual StatusCode sysExecute ()
virtual StatusCode preExecute ()
virtual StatusCode postExecute ()
virtual StatusCode sysFinalize ()
IDataProviderSvc * arcSvc () const
void putTES (DataObject *obj, std::string location) const
TYPE * getTES (std::string location) const
TYPE * getAES (std::string location, int index) const
std::vector< DataObject * > getAEScollection (std::string location) const
int getExecNum ()
std::string Location () const

Protected Member Functions

DayaBay::UserDataHeaderMakeHeaderObject ()
void InitializeHeader (DayaBay::HeaderObject *header)
TYPE * MakeHeader ()
TYPE * MakeHeader (std::vector< const DayaBay::IHeader * > &inputHeaders)
TYPE * MakeHeader (const DayaBay::IHeader *referenceHeader)
void AppendInputHeader (const DayaBay::HeaderObject *header) const

Protected Attributes

DayaBay::HeaderObjectm_headerObject
bool m_pullMode
std::string m_location

Private Attributes

std::string m_calibLocation
 Property CalibLocation - location in TES where the input CalibReadoutHeader is to be found.
std::vector< TimeStampm_lastTimestamp
 Property calibStatsLocation - location in TES where the output CalibStats is to be written.
std::vector< TimeStampm_lastMuonTimestamp
std::vector< TimeStampm_lastShowerTimestamp
std::vector< double > m_lastShowerEnergy
std::vector< TimeStampm_lastCalibShowerTimestamp
std::vector< double > m_lastCalibShowerEnergy
IDataProviderSvc * m_archiveSvc
double m_earlyTime
double m_lateTime
double m_peakEndTime
TH2F * m_mainPeak
TH2F * m_secondPeak
double m_nbins
double m_psd_nbins1
double m_psd_nbins2
TH2F * h_timepattern
TH1F * hpsd
TH1F * hpsd_local
TH1F * h_charge_column
int m_TotalBlockedTrigNum [DetectorId::kAll]
int m_TrigNum [DetectorId::kAll]
TimeStamp m_buffer_full_time_Sum [DetectorId::kAll]
TimeStamp m_BeginTime
std::string m_ChannelQualitySvcName
IChannelQualitySvcm_chanqualSvc
int nocqCount

Detailed Description

Definition at line 31 of file CalibStatsAlg.h.


Constructor & Destructor Documentation

CalibStatsAlg::CalibStatsAlg ( const std::string &  name,
ISvcLocator *  pSvcLocator 
)

Definition at line 57 of file CalibStatsAlg.cc.

  : DybAlgorithm<DayaBay::UserDataHeader>(name,pSvcLocator)
  , m_chanqualSvc(0)
{
  declareProperty("calibLocation", 
                  m_calibLocation=DayaBay::CalibReadoutHeaderLocation::Default,
                  "CalibReadoutHeader location in the TES.");
    declareProperty("earlyTime", m_earlyTime = -1650, "early time for charge count");
    declareProperty("lateTime", m_lateTime = -1250, "late time for charge count");
    declareProperty("peakEndTime", m_peakEndTime = -1480, "end time for peak charge count");
    declareProperty("nBINs", m_nbins = 800, "number of bins of the timing histogram (early to late)");
    declareProperty("PSD_nbin1",m_psd_nbins1 = 400, "number of bins for PSD1");
    declareProperty("PSD_nbin2",m_psd_nbins2 = 300, "number of bins for PSD2");
    declareProperty("ChannelQualitySvcName", m_ChannelQualitySvcName="DybChannelQualitySvc",
                    "Name of channel quality service");
}
CalibStatsAlg::~CalibStatsAlg ( ) [virtual]

Definition at line 74 of file CalibStatsAlg.cc.

{
}

Member Function Documentation

StatusCode CalibStatsAlg::initialize ( ) [virtual]

Definition at line 78 of file CalibStatsAlg.cc.

{
  StatusCode sc = this->GaudiAlgorithm::initialize();
  if(sc.isFailure()) return sc;
  m_lastTimestamp = std::vector<TimeStamp>(DetectorId::kAll, 
                                           TimeStamp::GetBOT());

  m_lastMuonTimestamp = std::vector<TimeStamp>(DetectorId::kAll, 
                                           TimeStamp::GetBOT());

  m_lastShowerTimestamp = std::vector<TimeStamp>(DetectorId::kAll, 
                                           TimeStamp::GetBOT());

  m_lastShowerEnergy = std::vector<double>(DetectorId::kAll, 
                                           0);
  // Added by Wenqiang, to tag shower by 2 inch calibPMT
  m_lastCalibShowerTimestamp = std::vector<TimeStamp>(DetectorId::kAll, 
                    TimeStamp::GetBOT());
  m_lastCalibShowerEnergy = std::vector<double>(DetectorId::kAll, 
                    0);
  //end of modification

  // Initialize livetime variables
  for(DetectorId::DetectorId_t detectorId=DetectorId::kAD1;
      detectorId!=DetectorId::kAll;
      detectorId=(DetectorId::DetectorId_t)(detectorId+1)){
    m_TotalBlockedTrigNum[detectorId]=0; 
    m_TrigNum[detectorId]=0;
    m_buffer_full_time_Sum[detectorId]=TimeStamp::GetBOT();
  }
  m_BeginTime=TimeStamp::GetBOT();

  // Get Archive event store
  sc = service("EventDataArchiveSvc",m_archiveSvc);
  if(sc.isFailure()){
    error() << "Failed to access the Archive event store" << endreq;
    return sc;
  }
  m_mainPeak = new TH2F("mainPeak", "", 13, 0, 13, 10, 0, 10);
  m_secondPeak = new TH2F("secondPeak", "", 13, 0, 13, 10, 0, 10);
        
  //added by XQ
  h_timepattern = new TH2F("timepattern","timepattern",24,1,25,8,1,9);
  hpsd = new TH1F("hpsd","hpsd",Int_t(m_nbins),m_earlyTime,m_lateTime);
  h_charge_column = new TH1F("h_charge_column","h_charge_column",24,1,25);
  hpsd_local = new TH1F("hpsd_local","hpsd_local",Int_t(m_nbins),m_earlyTime,m_lateTime);
  //end of modification


  // channel quality service
  m_chanqualSvc = svc<IChannelQualitySvc>(m_ChannelQualitySvcName,true);
  nocqCount=0;

  return sc;

}
StatusCode CalibStatsAlg::execute ( ) [virtual]

Get ReadoutHeader

Definition at line 135 of file CalibStatsAlg.cc.

{
  StatusCode sc;

  // Get CalibReadoutHeader. If nonexistent, skip event.
  if (!exist<CalibReadoutHeader>(m_calibLocation))
    return StatusCode::SUCCESS;
  
  CalibReadoutHeader* calibHeader =
    getTES<CalibReadoutHeader>(m_calibLocation);

  // Create output user data header
  UserDataHeader* calibStats = MakeHeaderObject();


  Context context = calibHeader->context();
  calibStats->setContext(context);
  calibStats->setEarliest(calibHeader->earliest());
  calibStats->setLatest(calibHeader->latest());

  
  DetectorId::DetectorId_t detectorId = calibHeader->context().GetDetId();

  // Get Calibrated Readout
  const CalibReadout* calibReadout = calibHeader->calibReadout();
  const CalibReadoutPmtCrate* pmtReadout = 0;
  if(!calibReadout) {
    //warning() << "No CalibReadout retrieved this cycle." << endreq;
  }else{
    if(detectorId!=DetectorId::kRPC && detectorId!=DetectorId::kUnknown){
      pmtReadout = dynamic_cast<const CalibReadoutPmtCrate*>(calibReadout);
      if(!pmtReadout) {
        error() << "Invalid PMT Crate readout!." << endreq;
        return StatusCode::FAILURE;
      }
    }
  }

  

  int triggerNumber=calibReadout->triggerNumber();
  int nHit=0;
  float nPESum=0;
  // Add a 2 inch PMT charge sum: ncalibPESum
  // by Wenqiang, wenqiang.gu@sjtu.edu.cn
  float nCalibPESum=0;
  // end of modification
  std::vector<float> calibPMTCharge (6,-1000);// initialize with a large negative number
  std::vector<int> calibPMTNHit (6,0);// 2 inch PMT nHit, initialize with a large negative number
  float nPEMedian=0;
  float nPERMS=0;
  float nPEMax=0;
  float nPulseSum=0;
  float nPulseMedian=0;
  float nPulseRMS=0;
  float tEarliest=0;
  float tLatest=0;
  float tMean=0;
  float tMedian=0;
  float tRMS=0;
  float earlyCharge=0;
  float lateCharge=0;
  float nominalCharge=0;
  float peakCharge=0;
  float nBadPMTs = 0;
  float nActivePMTs = 0;
  
  ServiceMode svcMode(context, 0);

  
  if(pmtReadout){
    
    // Get ChannelQualityService
    if (detectorId!=DetectorId::kIWS && detectorId!=DetectorId::kOWS){ // not used for WP for now
      const IChannelQuality * cq = m_chanqualSvc->channelQuality(svcMode);
      if (!cq){
        nocqCount ++;
        if (nocqCount <= NOCQ_ERRORS){
          warning() << "\tGot null channel quality object for " << context.AsString()
                    << ". Assuming all channels are good." << endreq;
        }
        if (nocqCount == NOCQ_ERRORS){
          warning() << "\tFurther warning for channel quality service will be suppressed" << endreq;
        }
        
        nBadPMTs = 0;
        nActivePMTs = 192;
      }else{
        // loop over channels and get number of bad PMTs
        IChannelQuality::ChannelSet_t chans = cq->channels();
        for (int i = 0; i < chans.size(); i++){
          DayaBay::FeeChannelId chid = chans[i];
          if (!cq->good(chid)){
            nBadPMTs += 1;
            //    std::cout <<  "\tGot a BAD one: " << chid.asString()
            //              <<  " board " << chid.board() 
            //              <<  " connector " << chid.connector() 
            //              <<  ", set = " << cq->hvRequested(chid)
            //              <<  ", get = " << cq->hvMeasured(chid)
            //              <<  std::endl;
          }else{
            nActivePMTs += 1;
          }
        }
      }
      //      std::cout << "----------- nBadPMTs = " << nBadPMTs << " / " << chans.size() << " ---------------" << std::endl;
    }
    
    
    
    static double nPEArr[MAX_STATS];
    static double nPulseArr[MAX_STATS];
    static double hitTimeArr[MAX_STATS];
    CalibReadoutPmtCrate::PmtChannelReadouts::const_iterator chanIter, 
      chanEnd = pmtReadout->channelReadout().end();
    for(chanIter=pmtReadout->channelReadout().begin(); chanIter != chanEnd; 
        chanIter++){
      if(nHit==MAX_STATS){
        error() << "More than " << MAX_STATS << " channels hit in this event."
                << " This is unphysical.  Bailing!" << endreq;
        return StatusCode::FAILURE;
      }
      const CalibReadoutPmtChannel& channel = *chanIter;

// Modified by Wenqiang, don't ignore 2 inch PMT, store charge sum
//      if(channel.pmtSensorId().isAD()){
//      AdPmtSensor pmtId(channel.pmtSensorId().fullPackedData());
//      if( pmtId.ring()==0 ){ continue; }  // Ignore calibration pmts
//      }
      if(channel.pmtSensorId().isAD()){
        AdPmtSensor pmtId(channel.pmtSensorId().fullPackedData());
        if( pmtId.ring()==0 ){ 
          // Get channel hit information
          std::vector<double> charges = channel.charge();
          std::vector<double> times = channel.time();

          // Calculate various charge sums
          float nCalibPE = 0;  // Charge sum of all hits for current PMT
          for(unsigned int i=0;i<times.size();i++){
            float hitCharge = charges[i];  // Current hit
            nCalibPE += hitCharge;  // All hits for this PMT
            //if(times[i]<m_earlyTime) earlyCharge += hitCharge; // all early hits
            //if(times[i]>m_lateTime) lateCharge += hitCharge; // all late hits
          }
          int col = pmtId.column(); //2inch PMT column 1 --> 6
          calibPMTCharge[col-1] = nCalibPE;
          calibPMTNHit[col-1] = times.size();
          nCalibPESum += nCalibPE;  // Add to total sum for event
          continue;
          
        }
      }
      // end of modification by Wenqiang

      // Get channel hit information
      std::vector<double> charges = channel.charge();
      std::vector<double> times = channel.time();

      // Calculate various charge sums
      float nPE = 0;  // Charge sum of all hits for current PMT
      for(unsigned int i=0;i<times.size();i++){
        float hitCharge = charges[i];  // Current hit
        nPE += hitCharge;  // All hits for this PMT
        if(times[i]<m_earlyTime) earlyCharge += hitCharge; // all early hits
        if(times[i]>m_lateTime) lateCharge += hitCharge; // all late hits
      }
      nPESum += nPE;  // Add to total sum for event

      // Nominal Charge Sum:
      //   This charge is our current 'best' estimate of the event
      //   charge.  This is equivalent to the charge sum used as input
      //   to the energy reconstruction.
      nominalCharge += channel.earliestCharge(m_earlyTime,m_lateTime);

      // Peak charge: sum of charge in the peak region only, as recommended by Soeren (Doc-8319)
      //   This charge is our current 'best' estimate of the event
      //   charge.  This is equivalent to the charge sum used as input
      //   to the energy reconstruction.
      peakCharge +=  channel.earliestCharge(m_earlyTime,m_peakEndTime);

      // General hit statistics
      float nPulses = channel.size();
      float hitTime = channel.earliestTime();
      if(nPE > nPEMax) nPEMax=nPE;
      nPEArr[nHit] = nPE;
      nPulseSum+=nPulses;
      nPulseArr[nHit] = nPulses;
      if(nHit==0){
        tEarliest = hitTime;
        tLatest = hitTime;
      }else{
        if(hitTime<tEarliest) tEarliest = hitTime;
        if(hitTime>tLatest) tLatest = hitTime;
      }
      tMean+=hitTime;
      hitTimeArr[nHit] = hitTime;
      nHit++;
    }
    if(nHit>0){
      tMean/=nHit;
      double nPEMean = nPESum/nHit;
      double nPulseMean = nPulseSum/nHit;
      for(int idx=0;idx<nHit;idx++){
        double deltaPE = nPEArr[idx] - nPEMean;
        nPERMS += (deltaPE*deltaPE);
        double deltaPulse = nPulseArr[idx] - nPulseMean;
        nPulseRMS += (deltaPulse*deltaPulse);
        double deltaTime = hitTimeArr[idx] - tMean;
        tRMS += (deltaTime*deltaTime);
      }
      nPERMS = sqrt(nPERMS/nHit);
      nPulseRMS = sqrt(nPulseRMS/nHit);
      tRMS = sqrt(tRMS/nHit);
      nPEMedian = TMath::Median(nHit, nPEArr);
      nPulseMedian = TMath::Median(nHit, nPulseArr);
      tMedian = TMath::Median(nHit, hitTimeArr);
    }
  }

  // =================== AD Flasher Discriminators
  double Quadrant = 0;
  double TimeRMS = 0;
  double MaxQ = 0;
  double MaxQ_2inchPMT = 0;
  int Column_2inchPMT = -1;
  double PeakRMS = 0;
  double Kurtosis = 0;
  int    FlasherPmtRing = 0;
  int    FlasherPmtColumn = 0;

  double quadrantQ1 = 0;
  double quadrantQ2 = 0;
  double quadrantQ3 = 0;
  double quadrantQ4 = 0;
  double ColumnKurtosis = 0;
  double RingKurtosis = 0;
  double MainPeakRMS = 0;
  double SecondPeakRMS = 0;

  // add by Xin 9/26/11 xqian@caltech.edu
  double time_PSD=0,time_PSD1=0;
  double time_PSD_local_RMS=0;
  double flasher_flag=0;
  double charge_sum_flasher_max=0;
  double Q1=0,Q2=0,Q3=0;
  int flasher_ring;
  int flasher_column;

   h_timepattern->Reset();
   hpsd->Reset();
   h_charge_column->Reset();
   hpsd_local->Reset();
   // end of modification


  if(pmtReadout && detectorId!=DetectorId::kIWS && detectorId!=DetectorId::kOWS){
    int    calibRing[MAX_STATS];
    int    calibColumn[MAX_STATS];
    double calibCharge[MAX_STATS];
    double calibHitCount[MAX_STATS];
    double calibTime[MAX_STATS];
    double channelMaxQ[nRing][nColumn];
    double pmtHit[nRing][nColumn];

    for(int iRing=0; iRing<nRing; iRing++){
      for(int iColumn=0; iColumn<nColumn; iColumn++){
        channelMaxQ[iRing][iColumn] = 0.;
        pmtHit[iRing][iColumn] = 0.;
      }
    }
    for(int i=0; i<MAX_STATS; i++) {
      calibRing[i] = 0;
      calibColumn[i] = 0;
      calibCharge[i] = 0.;
      calibHitCount[i] = 0.;
      calibTime[i] = 0.;
    }

    CalibReadoutPmtCrate::PmtChannelReadouts::const_iterator chanIter, 
      chanEnd = pmtReadout->channelReadout().end();
    int index = -1;
    int totalHitNum = 0; 
    for(chanIter=pmtReadout->channelReadout().begin(); chanIter != chanEnd; 
        chanIter++){
      const CalibReadoutPmtChannel& channel = *chanIter;
      AdPmtSensor pmtId(channel.pmtSensorId().fullPackedData());
      if(!pmtId.is8inch() && !pmtId.is2inch()) continue;
      int ring = -1;
      int clmn = -1;
      ring = pmtId.ring();
      clmn = pmtId.column();
      // ring: 0, column: 1-6 --> 2 inch PMT
      // ring: 1-8, column: 1-24 --> 8 inch PMT
      channelMaxQ[ring][clmn] = channel.maxCharge();
      debug() << "ring, clmn, channelMaxQ[ring][clmn]: " << ring << ",  " << clmn << ",  " << channelMaxQ[ring][clmn] << endreq;
      pmtHit[ring][clmn] = channel.size();
      totalHitNum += channel.size();

      // added by XQ
      double charge_sum_flasher = 0;
      double prev_charge = 0;
      double first_time = 0;
      double first_charge = 0;
      
      //end of modification

      for(unsigned int i=0;i<channel.size();i++){
        index++;
        
        if (first_time==0&&channel.time(i)>=m_earlyTime&&channel.time(i)<=m_lateTime) {
          first_time = channel.time(i);
          first_charge = channel.charge(i);
        }
        calibCharge[index] = channel.charge(i);
        calibTime[index] = channel.time(i);
        calibHitCount[index] = i;
        calibRing[index] = ring;
        calibColumn[index] = clmn;
        if (ring>0&&calibTime[index]>=m_earlyTime&&calibTime[index]<=m_lateTime)
          hpsd->Fill(calibTime[index]);
        

        if(!pmtId.is8inch()){
          if (channel.time(i)>=m_earlyTime&&channel.time(i)<=m_lateTime||1>0){
            if (channel.charge(i)>MaxQ_2inchPMT){
              MaxQ_2inchPMT = channel.charge(i);
              Column_2inchPMT = clmn;
            }
          }
          continue;   // Ignore calibration pmts
        }
        //added by XQ
        if (channel.charge(i)!=prev_charge&&channel.charge(i)>0.) charge_sum_flasher += channel.charge(i);
        prev_charge = channel.charge(i);
        //end of modification
      }
      
      //added by XQ
      if (index >=0&&clmn>0&&ring>0){
        h_timepattern->Fill(clmn,ring,first_time);
        h_charge_column->Fill(clmn,first_charge);
        if (charge_sum_flasher > charge_sum_flasher_max) 
          {
            charge_sum_flasher_max = charge_sum_flasher;
            flasher_ring = ring;
            flasher_column = clmn;
            //this might be the place to save more information about this so called flasher e.g. preADC, peakCycle etc.
          }
      }
      // end of modification
    }
    
    //================== Discriminators for 8 inch PMTs ========================
    int mainNhit = 0, secondNhit = 0;
    double effTime[192];
    int effNpmt = 0;
    for(int i=0; i<192; i++) effTime[i] = 0.;
  
    // calculate effective mean time and rms of all 8 inch PMTs
    for(int i=0; i<totalHitNum; i++) {
      if(calibTime[i]>m_earlyTime && calibTime[i]<m_lateTime && calibRing[i]>0){ // FIXME: select a time window
        AdPmtSensor pmtId(calibRing[i],calibColumn[i],
                          calibHeader->context().GetSite(),
                          calibHeader->context().GetDetId());
        if(!pmtId.is8inch()) continue;  // Protect from non-8inch signals
        int id = (calibRing[i]-1)*24 + calibColumn[i] - 1;
        if(effTime[id] > calibTime[i]) {
          effTime[id] = calibTime[i];
          debug() << "effTime[id]: "<< effTime[id] << endreq;
        }
      }
    }
    
    double meanTime = 0.;
    for(int i=0; i<192; i++) {
      if(effTime[i]!=0.) {  meanTime += effTime[i]; effNpmt++; }
    }
    meanTime = meanTime/effNpmt;
    TimeRMS = 0.;
    for(int i=0; i<192; i++) {
      if(effTime[i]!=0.) TimeRMS += (meanTime-effTime[i])*(meanTime-effTime[i]);
    }
    if(effNpmt>1) {  TimeRMS = sqrt(TimeRMS/(effNpmt-1)); }
  
    // =========================================================================
    // find largest charge of all 8 inch PMTs
    int maxqRing = -1, maxqCol = -1;
    double largest = 0.;
    double sum = 0.;
    // Chose the largest charge of all 8 inch PMTs in the readout window
    for(int iRing=1; iRing<nRing; iRing++){
      for(int iColumn=1; iColumn<nColumn; iColumn++){
        debug() << "iRing, iCol, channelMaxQ[ring][clmn]: " << iRing << ",  " << iColumn << ",  " << channelMaxQ[iRing][iColumn] << endreq;
        if(channelMaxQ[iRing][iColumn] > largest){
           largest = channelMaxQ[iRing][iColumn];
             maxqRing = iRing;
             maxqCol = iColumn;
          }
      }
    }
    debug() << "maxqRing, maxqCol, largest : " << maxqRing << ",  " << maxqCol << ", " << largest << endreq;
    // sum of the maxCharge of all 8 inch PMTs
    for(int iRing=1; iRing<nRing; iRing++){
      for(int iColumn=1; iColumn<nColumn; iColumn++){
        sum += channelMaxQ[iRing][iColumn];
      }
    }
    MaxQ = largest/sum;
    debug() << "MaxQ, sum: " << MaxQ << ",  " << sum << endreq;
  
    // =========================================================================
    // calculate the quadrant Q of 8 inch PMTs
    int col1, col2, col3, col4, col5, col6;
    col3 = maxqCol;
    col1 = (col3-2)>0 ? (col3-2) : (col3-2) + 24;
    col2 = (col3-1)>0 ? (col3-1) : (col3-1) + 24;
    col4 = (col3+1)>24 ? (col3+1) - 24 : (col3+1);
    col5 = (col3+2)>24 ? (col3+2) - 24 : (col3+2);
    col6 = (col3+3)>24 ? (col3+3) - 24 : (col3+3);
    quadrantQ1 = 0.;
    for(int iRing=1; iRing<=8; iRing++) {
      quadrantQ1 += channelMaxQ[iRing][col1];
      quadrantQ1 += channelMaxQ[iRing][col2];
      quadrantQ1 += channelMaxQ[iRing][col3];
      quadrantQ1 += channelMaxQ[iRing][col4];
      quadrantQ1 += channelMaxQ[iRing][col5];
      quadrantQ1 += channelMaxQ[iRing][col6];
    }
  
    col3 = (maxqCol+12)>24 ? (maxqCol+12)-24 : (maxqCol+12); 
    col1 = (col3-2)>0 ? (col3-2) : (col3-2) + 24;
    col2 = (col3-1)>0 ? (col3-1) : (col3-1) + 24;
    col4 = (col3+1)>24 ? (col3+1) - 24 : (col3+1);
    col5 = (col3+2)>24 ? (col3+2) - 24 : (col3+2);
    col6 = (col3+3)>24 ? (col3+3) - 24 : (col3+3);
    quadrantQ3 = 0.;
    for(int iRing=1; iRing<=8; iRing++) {
      quadrantQ3 += channelMaxQ[iRing][col1];
      quadrantQ3 += channelMaxQ[iRing][col2];
      quadrantQ3 += channelMaxQ[iRing][col3];
      quadrantQ3 += channelMaxQ[iRing][col4];
      quadrantQ3 += channelMaxQ[iRing][col5];
      quadrantQ3 += channelMaxQ[iRing][col6];
    }
  
    col3 = (maxqCol+6)>24 ? (maxqCol+6)-24 : (maxqCol+6); 
    col1 = (col3-2)>0 ? (col3-2) : (col3-2) + 24;
    col2 = (col3-1)>0 ? (col3-1) : (col3-1) + 24;
    col4 = (col3+1)>24 ? (col3+1) - 24 : (col3+1);
    col5 = (col3+2)>24 ? (col3+2) - 24 : (col3+2);
    col6 = (col3+3)>24 ? (col3+3) - 24 : (col3+3);
    quadrantQ2 = 0.;
    for(int iRing=1; iRing<=8; iRing++) {
      quadrantQ2 += channelMaxQ[iRing][col1];
      quadrantQ2 += channelMaxQ[iRing][col2];
      quadrantQ2 += channelMaxQ[iRing][col3];
      quadrantQ2 += channelMaxQ[iRing][col4];
      quadrantQ2 += channelMaxQ[iRing][col5];
      quadrantQ2 += channelMaxQ[iRing][col6];
    }
  
    col3 = (maxqCol+18)>24 ? (maxqCol+18)-24 : (maxqCol+18); 
    col1 = (col3-2)>0 ? (col3-2) : (col3-2) + 24;
    col2 = (col3-1)>0 ? (col3-1) : (col3-1) + 24;
    col4 = (col3+1)>24 ? (col3+1) - 24 : (col3+1);
    col5 = (col3+2)>24 ? (col3+2) - 24 : (col3+2);
    col6 = (col3+3)>24 ? (col3+3) - 24 : (col3+3);
    quadrantQ4 = 0.;
    for(int iRing=1; iRing<=8; iRing++) {
      quadrantQ4 += channelMaxQ[iRing][col1];
      quadrantQ4 += channelMaxQ[iRing][col2];
      quadrantQ4 += channelMaxQ[iRing][col3];
      quadrantQ4 += channelMaxQ[iRing][col4];
      quadrantQ4 += channelMaxQ[iRing][col5];
      quadrantQ4 += channelMaxQ[iRing][col6];
    }
  
    if((quadrantQ2+quadrantQ4)!=0.) {
      Quadrant = quadrantQ3/(quadrantQ2+quadrantQ4);
    }
    else { // Some retriggers after muon has many hits with zero charges
      Quadrant = 10.;
    }
    //=========================================================================
    // calculate the kurtosis of two dimensional charge pattern
    double sum1x,sum2x, sum1y, sum2y;
    sum1x = sum2x = sum1y = sum2y = 0;
    double qsum1, qsum2;
    qsum1 = qsum2 = 0;
    double average1x, average2x, average1y, average2y;
    average1x = average2x = average1y = average2y = 0;
    double sd1x, sd2x, sd1y, sd2y;
    sd1x = sd2x = sd1y = sd2y = 0;
    double kurtosis1x, kurtosis2x, kurtosis1y, kurtosis2y;
    kurtosis1x = kurtosis2x = kurtosis1y = kurtosis2y = 0;
    // calculate the average
    for(int iRing=1; iRing<nRing; iRing++){
      for(int iColumn=1; iColumn<nColumn; iColumn++){
        int iCol = iColumn - maxqCol + 6;
        if(iCol>24) iCol = iCol - 24;
        if(iCol<1) iCol = iCol + 24;
        if(iCol<=12){
          m_mainPeak->SetBinContent(iCol+1, iRing+1, channelMaxQ[iRing][iColumn]);
          mainNhit++;
          qsum1 = qsum1 + channelMaxQ[iRing][iColumn];
          sum1x = sum1x + iCol*channelMaxQ[iRing][iColumn];
          sum1y = sum1y + iRing*channelMaxQ[iRing][iColumn];
        }
        if(iCol>12){
          m_secondPeak->SetBinContent(iCol+1-12, iRing+1, channelMaxQ[iRing][iColumn]);
          secondNhit++;
          qsum2 = qsum2 + channelMaxQ[iRing][iColumn];
          sum2x = sum2x + (iCol-12)*channelMaxQ[iRing][iColumn];
          sum2y = sum2y + iRing*channelMaxQ[iRing][iColumn];
        }
      }
    }
        average1x = sum1x/qsum1;
        average1y = sum1y/qsum1;
        average2x = sum2x/qsum2;
        average2y = sum2y/qsum2;
  
    // calculate the RMS 
    for(int iRing=1; iRing<nRing; iRing++){
      for(int iColumn=1; iColumn<nColumn; iColumn++){
        int iCol = iColumn - maxqCol + 6;
        if(iCol>24) iCol = iCol - 24;
        if(iCol<1) iCol = iCol + 24;
        if(iCol<=12 && channelMaxQ[iRing][iColumn]>0){
          sd1x = sd1x + channelMaxQ[iRing][iColumn]*pow(iCol-average1x, 2);
          sd1y = sd1y + channelMaxQ[iRing][iColumn]*pow(iRing-average1y, 2);
        }
        if(iCol>12 && channelMaxQ[iRing][iColumn]>0){
          sd2x = sd2x + channelMaxQ[iRing][iColumn]*pow(iCol-12-average2x, 2);
          sd2y = sd2y + channelMaxQ[iRing][iColumn]*pow(iRing-average2y, 2);
        }
      }
    }
        sd1x = sqrt(sd1x/qsum1);
        sd1y = sqrt(sd1y/qsum1);
        sd2x = sqrt(sd2x/qsum2);
        sd2y = sqrt(sd2y/qsum2);
  
    // calculate the kurtosis
    for(int iRing=1; iRing<nRing; iRing++){
      for(int iColumn=1; iColumn<nColumn; iColumn++){
        int iCol = iColumn - maxqCol + 6;
        if(iCol>24) iCol = iCol - 24;
        if(iCol<1) iCol = iCol + 24;
        if(iCol<=12 && channelMaxQ[iRing][iColumn]>0){
          kurtosis1x = kurtosis1x + channelMaxQ[iRing][iColumn]*pow(iCol-average1x, 4);
          kurtosis1y = kurtosis1y + channelMaxQ[iRing][iColumn]*pow(iRing-average1y, 4);
        }
        if(iCol>12 && channelMaxQ[iRing][iColumn]>0){
          kurtosis2x = kurtosis2x + channelMaxQ[iRing][iColumn]*pow(iCol-12-average2x, 4);
          kurtosis2y = kurtosis2y + channelMaxQ[iRing][iColumn]*pow(iRing-average2y, 4);
        }
      }
    }
    kurtosis1x = kurtosis1x/pow(sd1x, 4)/qsum1 - 3;
    kurtosis1y = kurtosis1y/pow(sd1y, 4)/qsum1 - 3;
    kurtosis2x = kurtosis2x/pow(sd2x, 4)/qsum2 - 3;
    kurtosis2y = kurtosis2y/pow(sd2y, 4)/qsum2 - 3;
  
    MainPeakRMS = sqrt(pow(m_mainPeak->GetRMS(1),2) + pow(m_mainPeak->GetRMS(2),2));
    SecondPeakRMS = sqrt(pow(m_secondPeak->GetRMS(1),2) + pow(m_secondPeak->GetRMS(2),2));
    PeakRMS = MainPeakRMS + SecondPeakRMS;
    ColumnKurtosis = kurtosis1x;
    RingKurtosis = kurtosis1y;
    Kurtosis = kurtosis1x + kurtosis1y;
    FlasherPmtRing = maxqRing;
    FlasherPmtColumn = maxqCol;
  }
  
  //added by XQ
  //calculate the local trms for four different methods 
  for(int k = flasher_column-2;k<=flasher_column+2;k++){
    for (int z = flasher_ring-2;z<=flasher_ring+2;z++){
      int kk=k;
      if (k>24) kk=k-24;
      if (k<1) kk=k+24;
      if (kk>=1&&kk<=24&&z>=1&&z<=8){
        if (h_timepattern->GetBinContent(kk,z)>=m_earlyTime&&h_timepattern->GetBinContent(kk,z)<=m_lateTime)
          hpsd_local->Fill(h_timepattern->GetBinContent(kk,z));
      }
    }
  }

  for(Int_t k=1;k!=25;k++){
    double tempcharge = h_charge_column->GetBinContent(k);
    double dis = getdis((double)k,(double)flasher_column);
    if (dis<=3) Q1 += tempcharge;
    if (dis>3&&dis<=9) Q2 += tempcharge;
    if (dis>9) Q3 += tempcharge;
  }
  
  
  time_PSD_local_RMS = hpsd_local->GetRMS();
  if (hpsd->GetSum()>0){
    time_PSD = hpsd->Integral(1,Int_t(m_psd_nbins1))/hpsd->GetSum();
    time_PSD1 = hpsd->Integral(1,Int_t(m_psd_nbins2))/hpsd->GetSum();
    if (Q2>0.){
      flasher_flag = nPEMax/nPESum + Q3/Q2*2 + (1-time_PSD) + (1-time_PSD1) + time_PSD_local_RMS/100.;
    }else{
      flasher_flag = nPEMax/nPESum + 1. + (1-time_PSD) + (1-time_PSD1) + time_PSD_local_RMS/100.;
    }
  }
  //std::cout << tRMS_local << "\t" << h_time->GetSum() << std::endl;
  //end of modification

  // Fill the statistics
  calibStats->set("triggerNumber",triggerNumber);
  calibStats->set("nHit",nHit);
  calibStats->set("nPESum",nPESum);
  //Added by Wenqiang
  calibStats->set("nCalibPESum",nCalibPESum);
  //
  calibStats->set("calibPMTCharge",calibPMTCharge);
  calibStats->set("calibPMTNHit",calibPMTNHit);
  calibStats->set("nPEMedian",nPEMedian);
  calibStats->set("nPERMS",nPERMS);
  calibStats->set("nPEMax",nPEMax);
  calibStats->set("nPulseSum",nPulseSum);
  calibStats->set("nPulseMedian",nPulseMedian);
  calibStats->set("nPulseRMS",nPulseRMS);
  calibStats->set("tEarliest",tEarliest);
  calibStats->set("tLatest",tLatest);
  calibStats->set("tMean",tMean);
  calibStats->set("tMedian",tMedian);
  calibStats->set("tRMS",tRMS);
  calibStats->set("EarlyCharge",earlyCharge);
  calibStats->set("LateCharge",lateCharge);
  calibStats->set("NominalCharge",nominalCharge);
  calibStats->set("PeakCharge",peakCharge);

  // number of bad/active PMTs from ChannelQualityService
  calibStats->set("nActivePMTs",nActivePMTs); 
  calibStats->set("nBadPMTs",nBadPMTs); 

  
  //Fill the PMT Flasher discriminators
  calibStats->set("maxqRing", FlasherPmtRing);
  calibStats->set("maxqCol", FlasherPmtColumn);
  calibStats->set("MiddleTimeRMS", (float)TimeRMS);
  calibStats->set("MaxQ", (float)MaxQ);
  calibStats->set("MaxQ_2inchPMT", (float)MaxQ_2inchPMT);
  calibStats->set("Column_2inchPMT", Column_2inchPMT);

  calibStats->set("Quadrant", (float)Quadrant);
  calibStats->set("PeakRMS", (float)PeakRMS);
  calibStats->set("Kurtosis", (float)Kurtosis);

  calibStats->set("QuadrantQ1", (float)quadrantQ1);
  calibStats->set("QuadrantQ2", (float)quadrantQ2);
  calibStats->set("QuadrantQ3", (float)quadrantQ3);
  calibStats->set("QuadrantQ4", (float)quadrantQ4);
  calibStats->set("RingKurtosis", (float)RingKurtosis);
  calibStats->set("ColumnKurtosis", (float)ColumnKurtosis);
  calibStats->set("MainPeakRMS", (float)MainPeakRMS);
  calibStats->set("SecondPeakRMS", (float)SecondPeakRMS);

  // add some additional flasher variable XQ

  calibStats->set("charge_sum_flasher_max",(float)charge_sum_flasher_max);
  calibStats->set("time_PSD",(float)time_PSD);
  calibStats->set("time_PSD1",(float)time_PSD1);
  calibStats->set("flasher_flag",(float)flasher_flag);
  calibStats->set("flasher_ring",flasher_ring);
  calibStats->set("flasher_column",flasher_column);
  calibStats->set("time_PSD_local_RMS",(float)time_PSD_local_RMS);
  calibStats->set("Q1",(float)Q1);
  calibStats->set("Q2",(float)Q2);
  calibStats->set("Q3",(float)Q3);
  
  //end of modification

  //info() << "dt loop" << endreq;
  // Time since previous triggers info
  for(DetectorId::DetectorId_t lastDetectorId=DetectorId::kAD1;
      lastDetectorId!=DetectorId::kAll;
      lastDetectorId=(DetectorId::DetectorId_t)(lastDetectorId+1)){
    float dtLast = -1.0*Units::ms;
    //info() << "   detId " << (int)lastDetectorId << endreq; 
    if(m_lastTimestamp[lastDetectorId]!=TimeStamp::GetBOT()){
      TimeStamp dtLastTS(calibHeader->timeStamp());
      dtLastTS.Subtract(m_lastTimestamp[lastDetectorId]);
      dtLast = dtLastTS.GetSeconds()*Units::second;
    }
    calibStats->set(
     TString::Format("dtLast%s_ms",DetectorId::AsString(lastDetectorId)).Data(),
     float(dtLast / Units::ms));
  }

  //info() << "-----------"<<endreq;
  //info() << "-----------"<<endreq;
  // Addition from Bryce, 10-31: Time since AD muon and shower events, and energy of last shower
  for(DetectorId::DetectorId_t lastDetectorId=DetectorId::kAD1;
      lastDetectorId<(DetectorId::kIWS);
      lastDetectorId=(DetectorId::DetectorId_t)(lastDetectorId+1)){
    float dtLastMuon = -1.0*Units::ms;
    float dtLastShower = -1.0*Units::ms;
    float ELastShower = 0;
    // Two variables added by Wenqiang
    float dtLastCalibShower = -1.0*Units::ms;
    float ELastCalibShower = 0;
    // end of modification
    if(m_lastMuonTimestamp[lastDetectorId]!=TimeStamp::GetBOT()){
      //See how long since last muon
      TimeStamp dtLastMuonTS(calibHeader->timeStamp());
      dtLastMuonTS.Subtract(m_lastMuonTimestamp[lastDetectorId]);
      dtLastMuon = dtLastMuonTS.GetSeconds()*Units::second;
    }
    if(m_lastShowerTimestamp[lastDetectorId]!=TimeStamp::GetBOT()){
      //See how long since last shower
      TimeStamp dtLastShowerTS(calibHeader->timeStamp());
      dtLastShowerTS.Subtract(m_lastShowerTimestamp[lastDetectorId]);
      dtLastShower = dtLastShowerTS.GetSeconds()*Units::second;
      ELastShower = m_lastShowerEnergy[lastDetectorId];
    }
    // Added by Wenqiang
    if(m_lastCalibShowerTimestamp[lastDetectorId]!=TimeStamp::GetBOT()){
      TimeStamp dtLastCalibShowerTS(calibHeader->timeStamp());
      dtLastCalibShowerTS.Subtract(m_lastCalibShowerTimestamp[lastDetectorId]);
      dtLastCalibShower = dtLastCalibShowerTS.GetSeconds()*Units::second;
      ELastCalibShower = m_lastCalibShowerEnergy[lastDetectorId];
      // End of modification
    }

//    else{
//      m_lastMuonTimestamp[lastDetectorId] = calibHeader->timeStamp();
//      m_lastShowerTimestamp[lastDetectorId] = calibHeader->timeStamp();
//      m_lastShowerEnergy[lastDetectorId] = ELastShower;
//    }      
    
//    info() << "      " << "Detector type: "<<TString::Format("%s",DetectorId::AsString(detectorId)).Data()<<endreq;
//    info() << "      " << "dtLast_ADMuon_ms " <<dtLastMuon/Units::ms<< endreq; 
//    info() << "      " << "dtLast_ADShower_ms " <<dtLastShower/Units::ms<< endreq; 
//    info() << "      " << "ELast_ADShower_pe " <<ELastShower<< endreq; 
    if(detectorId>=(DetectorId::kIWS)){
      //info() << "Wrote from"<<TString::Format("%s",DetectorId::AsString(detectorId)).Data()<<endreq;
      calibStats->set("dtLast_ADMuon_ms",0);
      calibStats->set("dtLast_ADShower_ms",0);
      calibStats->set("ELast_ADShower_pe",0);
      // Added by Wenqiang
      calibStats->set("dtLast_ADCalibShower_ms",0);
      calibStats->set("ELast_ADCalibShower_pe",0);
      // end of modification
    }
    if(lastDetectorId == detectorId){
      //info() << "Wrote from"<<TString::Format("%s",DetectorId::AsString(lastDetectorId)).Data()<<endreq;
      calibStats->set("dtLast_ADMuon_ms",float(dtLastMuon / Units::ms));
      calibStats->set("dtLast_ADShower_ms",float(dtLastShower / Units::ms));
      calibStats->set("ELast_ADShower_pe",float(ELastShower));
      // Added by Wenqiang
      calibStats->set("dtLast_ADCalibShower_ms",
                    float(dtLastCalibShower / Units::ms));
      calibStats->set("ELast_ADCalibShower_pe",float(ELastCalibShower));
      // end of modification
    }
    //info() << "-----------"<<endreq;
  }
  //10-31: End of Bryce's edit.
    

  // Time until next triggers info (placeholder)
  for(DetectorId::DetectorId_t nextDetectorId=DetectorId::kAD1;
      nextDetectorId!=DetectorId::kAll;
      nextDetectorId=(DetectorId::DetectorId_t)(nextDetectorId+1)){
    calibStats->set(
     TString::Format("dtNext%s_ms",DetectorId::AsString(nextDetectorId)).Data(),
     float(-1.0));    
  }

  // Try to fill dtNext information for previous triggers (in AES)
  if(m_archiveSvc 
     && exist<DybArchiveList>(m_archiveSvc,m_location)){
    DybArchiveList* arcList = 0;
    arcList = get<DybArchiveList*>(m_archiveSvc,m_location);
    if(!arcList){
      error() << "Failed to get archive of CalibStats" << endreq;
      return StatusCode::FAILURE;
    }
    for(unsigned int headerIdx=0; headerIdx<arcList->size(); headerIdx++){
      UserDataHeader* archiveCalibStats 
        = dynamic_cast<UserDataHeader*>(arcList->dataObject(headerIdx));
      if(!archiveCalibStats){
        error() << "Invalid CalibStats in AES!" << endreq;
        return StatusCode::FAILURE;
      }
      TimeStamp dtNextTS(calibHeader->timeStamp());
      dtNextTS.Subtract(archiveCalibStats->timeStamp());
      //info() << TString::Format("dtNext%s_ms",DetectorId::AsString(detectorId)).Data() << " = " << float(dtNextTS.GetSeconds()*Units::second/Units::ms)
      //             << endreq;
      archiveCalibStats->set(
         TString::Format("dtNext%s_ms",DetectorId::AsString(detectorId)).Data(),
         float(dtNextTS.GetSeconds()*Units::second/Units::ms));
      if(archiveCalibStats->context().GetDetId()==detectorId){
        //Same detector as current CalibStats, should not need to look
        //further in archive.
        break;
      }
    }
  }


  //10-31 Bryce edit
  //keep track of last AD Muon and Shower Muon
  if(nominalCharge > 3200)
    m_lastMuonTimestamp[detectorId] = calibHeader->timeStamp();
  if(nominalCharge > 160000){
    m_lastShowerTimestamp[detectorId] = calibHeader->timeStamp();
    m_lastShowerEnergy[detectorId] = nominalCharge;
  }
   // Added by wenqiang, keep track of 2 inch calibPMT tagged shower
  if(nCalibPESum > 200){
    m_lastCalibShowerTimestamp[detectorId] = calibHeader->timeStamp();
    m_lastCalibShowerEnergy[detectorId] = nCalibPESum;
  }
  // End of modification
 
  // Added by BeiZhen for runtime and livetime estimate: ===================
  //  Set placeholder variables for CalibStats data 
  int buffer_full_flag = -1;
  int nBlockedTriggers = -1;
  float integralLiveTime_block_trigger_ms = -1;
  float integralLiveTime_buffer_full_ms = -1;
  float integralRunTime_ms = -1;
  // place holder for variable saving total number of channels having readouts
  int nTotalChannels = -1;

  // First, calculate integral runtime
  if( m_BeginTime==TimeStamp::GetBOT() ){
    // First trigger
    m_BeginTime = calibReadout->triggerTime();
    integralRunTime_ms = 0;
  }else{
    TimeStamp integralRunTS = calibReadout->triggerTime();
    integralRunTS.Subtract(m_BeginTime);
    integralRunTime_ms = integralRunTS.GetSeconds()*Units::second/Units::ms;
  }

  // Next, calculate integral livetime, and related quantities
  if(detectorId!=DetectorId::kRPC 
     && detectorId!=DetectorId::kAll 
     && detectorId!=DetectorId::kUnknown){
    m_TrigNum[detectorId]++;

    const DayaBay::DaqLtbFrame* firstTriggerFrame = 0;
    ReadoutHeader* readoutHdr = 
      get<ReadoutHeader>(ReadoutHeader::defaultLocation());
    if(readoutHdr) {
      const DaqCrate* daqCrate = readoutHdr->daqCrate();
      if(daqCrate){
        const DaqPmtCrate* pmtCrate = daqCrate->asPmtCrate();
        if(pmtCrate){
          const DaqLtb& pmtLtb = pmtCrate->localTriggerBoard();
          const DayaBay::DaqLtb::LtbFramePtrList& pmtLtbFramesPtrs 
            = pmtLtb.frames();
          firstTriggerFrame = pmtLtbFramesPtrs.front();
          nTotalChannels = pmtCrate->channelReadouts().size();
        }
      }
    }

    if(!firstTriggerFrame){
      static int warningCount=0;
      if(warningCount<10){
        warning() << "Failed to get raw LTB data;"
                  << " livetime variables will not be correct"<<endreq;
        warningCount++;
        if(warningCount==10){
          warning() << "Further 'readout header' warnings will suppressed."
                    <<endreq;
        }
      }
    }else{
      // Found LTB data, calculate livetime
      
      nBlockedTriggers = firstTriggerFrame->blockedTriggerCount();
      // total # of blocked Trigger for each detector
      m_TotalBlockedTrigNum[detectorId] += nBlockedTriggers; 
      buffer_full_flag = (firstTriggerFrame->feeBufferFull() 
                          || firstTriggerFrame->ltbBufferFull());
      if( buffer_full_flag && (nBlockedTriggers > 0) ){
        // Missing some triggers; assume dead time is approximately time
        // since last trigger in this same detector.
        if(m_lastTimestamp[detectorId]!=TimeStamp::GetBOT()){
          TimeStamp dtLastTS(calibHeader->timeStamp());
          dtLastTS.Subtract(m_lastTimestamp[detectorId]);
          m_buffer_full_time_Sum[detectorId].Add(dtLastTS);
        }
      }
        
      double blockedRatio = (m_TotalBlockedTrigNum[detectorId]
                             /(m_TotalBlockedTrigNum[detectorId]
                               +m_TrigNum[detectorId]));
      integralLiveTime_block_trigger_ms = integralRunTime_ms * (1 - blockedRatio);
      integralLiveTime_buffer_full_ms =
        (integralRunTime_ms
         - (m_buffer_full_time_Sum[detectorId].GetSeconds()
            *Units::second/Units::ms));
    }
  }
  // Set livetime related variables
  calibStats->set("blocked_trigger",nBlockedTriggers);
  calibStats->set("integralLiveTime_blocked_trigger_ms",
                  integralLiveTime_block_trigger_ms);
  calibStats->set("integralLiveTime_buffer_full_ms",
                  integralLiveTime_buffer_full_ms);
  calibStats->set("integralRunTime_ms",integralRunTime_ms);
  calibStats->set("buffer_full_flag",buffer_full_flag);
  // end runtime/livetime estimate
  calibStats->set("nTotalChannels",nTotalChannels);

  //Just keep track of last event
  m_lastTimestamp[detectorId] = calibHeader->timeStamp();
  
  return StatusCode::SUCCESS;
}
StatusCode CalibStatsAlg::finalize ( ) [virtual]

Definition at line 1056 of file CalibStatsAlg.cc.

{

  //Status report for null Channel Quality Service
  if (nocqCount > 0)
    warning() << "******* Got null Channel Quality objects for "<< nocqCount << " times  *********" << endreq;



  if(m_archiveSvc){
    m_archiveSvc->release();
    m_archiveSvc=0;
  }
  delete m_mainPeak;
  delete m_secondPeak;
  delete h_timepattern;
  delete hpsd;
  delete h_charge_column;
  delete hpsd_local;
  return this->GaudiAlgorithm::finalize();
}
DayaBay::UserDataHeader * DybAlgorithm< DayaBay::UserDataHeader >::GetCurrentHeaderObject ( ) const [inherited]
virtual StatusCode DybAlgorithm< DayaBay::UserDataHeader >::sysInitialize ( ) [virtual, inherited]

Reimplemented from DybBaseAlg.

virtual StatusCode DybAlgorithm< DayaBay::UserDataHeader >::sysExecute ( ) [virtual, inherited]

Reimplemented from DybBaseAlg.

virtual StatusCode DybAlgorithm< DayaBay::UserDataHeader >::preExecute ( ) [virtual, inherited]

Reimplemented from DybBaseAlg.

virtual StatusCode DybAlgorithm< DayaBay::UserDataHeader >::postExecute ( ) [virtual, inherited]

Reimplemented from DybBaseAlg.

virtual StatusCode DybAlgorithm< DayaBay::UserDataHeader >::sysFinalize ( ) [virtual, inherited]

Reimplemented from DybBaseAlg.

IDataProviderSvc * DybAlgorithm< DayaBay::UserDataHeader >::arcSvc ( ) const [inherited]

Reimplemented from DybBaseAlg.

void DybAlgorithm< DayaBay::UserDataHeader >::putTES ( DataObject *  obj,
std::string  location 
) const [inherited]

Reimplemented from DybBaseAlg.

TYPE * DybAlgorithm< DayaBay::UserDataHeader >::getTES ( std::string  location) const [inherited]

Reimplemented from DybBaseAlg.

TYPE * DybAlgorithm< DayaBay::UserDataHeader >::getAES ( std::string  location,
int  index 
) const [inherited]

Reimplemented from DybBaseAlg.

std::vector< DataObject * > DybAlgorithm< DayaBay::UserDataHeader >::getAEScollection ( std::string  location) const [inherited]

Reimplemented from DybBaseAlg.

int DybAlgorithm< DayaBay::UserDataHeader >::getExecNum ( ) [inherited]

Reimplemented from DybBaseAlg.

std::string DybAlgorithm< DayaBay::UserDataHeader >::Location ( ) const [inherited]

Reimplemented from DybBaseAlg.

DayaBay::UserDataHeader * DybAlgorithm< DayaBay::UserDataHeader >::MakeHeaderObject ( ) [protected, inherited]
void DybAlgorithm< DayaBay::UserDataHeader >::InitializeHeader ( DayaBay::HeaderObject header) [protected, inherited]

Reimplemented from DybBaseAlg.

TYPE * DybAlgorithm< DayaBay::UserDataHeader >::MakeHeader ( ) [protected, inherited]

Reimplemented from DybBaseAlg.

TYPE * DybAlgorithm< DayaBay::UserDataHeader >::MakeHeader ( std::vector< const DayaBay::IHeader * > &  inputHeaders) [protected, inherited]

Reimplemented from DybBaseAlg.

TYPE * DybAlgorithm< DayaBay::UserDataHeader >::MakeHeader ( const DayaBay::IHeader referenceHeader) [protected, inherited]

Reimplemented from DybBaseAlg.

void DybAlgorithm< DayaBay::UserDataHeader >::AppendInputHeader ( const DayaBay::HeaderObject header) const [protected, inherited]

Reimplemented from DybBaseAlg.


Member Data Documentation

std::string CalibStatsAlg::m_calibLocation [private]

Property CalibLocation - location in TES where the input CalibReadoutHeader is to be found.

Default is DayaBay::CalibReadoutHeaderLocation::Default

Definition at line 47 of file CalibStatsAlg.h.

std::vector<TimeStamp> CalibStatsAlg::m_lastTimestamp [private]

Property calibStatsLocation - location in TES where the output CalibStats is to be written.

Default is "/Event/Data/CalibStats"

Definition at line 56 of file CalibStatsAlg.h.

Definition at line 57 of file CalibStatsAlg.h.

Definition at line 58 of file CalibStatsAlg.h.

std::vector<double> CalibStatsAlg::m_lastShowerEnergy [private]

Definition at line 59 of file CalibStatsAlg.h.

Definition at line 61 of file CalibStatsAlg.h.

std::vector<double> CalibStatsAlg::m_lastCalibShowerEnergy [private]

Definition at line 62 of file CalibStatsAlg.h.

IDataProviderSvc* CalibStatsAlg::m_archiveSvc [private]

Definition at line 65 of file CalibStatsAlg.h.

double CalibStatsAlg::m_earlyTime [private]

Definition at line 66 of file CalibStatsAlg.h.

double CalibStatsAlg::m_lateTime [private]

Definition at line 67 of file CalibStatsAlg.h.

double CalibStatsAlg::m_peakEndTime [private]

Definition at line 68 of file CalibStatsAlg.h.

TH2F* CalibStatsAlg::m_mainPeak [private]

Definition at line 69 of file CalibStatsAlg.h.

TH2F* CalibStatsAlg::m_secondPeak [private]

Definition at line 70 of file CalibStatsAlg.h.

double CalibStatsAlg::m_nbins [private]

Definition at line 74 of file CalibStatsAlg.h.

double CalibStatsAlg::m_psd_nbins1 [private]

Definition at line 75 of file CalibStatsAlg.h.

double CalibStatsAlg::m_psd_nbins2 [private]

Definition at line 76 of file CalibStatsAlg.h.

Definition at line 78 of file CalibStatsAlg.h.

TH1F* CalibStatsAlg::hpsd [private]

Definition at line 79 of file CalibStatsAlg.h.

TH1F* CalibStatsAlg::hpsd_local [private]

Definition at line 80 of file CalibStatsAlg.h.

Definition at line 81 of file CalibStatsAlg.h.

Definition at line 86 of file CalibStatsAlg.h.

Definition at line 87 of file CalibStatsAlg.h.

Definition at line 88 of file CalibStatsAlg.h.

Definition at line 89 of file CalibStatsAlg.h.

Definition at line 94 of file CalibStatsAlg.h.

Definition at line 96 of file CalibStatsAlg.h.

int CalibStatsAlg::nocqCount [private]

Definition at line 97 of file CalibStatsAlg.h.

Reimplemented from DybBaseAlg.

bool DybAlgorithm< DayaBay::UserDataHeader >::m_pullMode [protected, inherited]

Reimplemented from DybBaseAlg.

std::string DybAlgorithm< DayaBay::UserDataHeader >::m_location [protected, inherited]

Reimplemented from DybBaseAlg.


The documentation for this class was generated from the following files:
| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |

Generated on Fri May 16 2014 10:13:52 for CalibStats by doxygen 1.7.4