/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 | Private Types | Private Member Functions | Private Attributes
ReconDataHistogram Class Reference

#include <ReconDataHistogram.h>

Collaboration diagram for ReconDataHistogram:
Collaboration graph
[legend]

List of all members.

Public Member Functions

 ReconDataHistogram (const std::string &name, ISvcLocator *pSvcLocator)
virtual ~ReconDataHistogram ()
virtual StatusCode initialize ()
virtual StatusCode finalize ()
virtual StatusCode execute ()

Private Types

enum  Fig_t {
  kEnergyVsR, kEnergyVsZ, kRateVsR, kRateVsZ,
  kEnergy, kZVsR, kYVsX
}

Private Member Functions

TH1 * getOrMakeHist (int run, const DayaBay::Detector &detector, int histIndex)
std::string getPath (int run, const DayaBay::Detector &detector, const char *histName)

Private Attributes

IStatisticsSvc * m_statsSvc
TimeStampm_firstTriggerTime
TimeStamp m_lastTriggerTime
std::map< int, TH1 ** > m_hist
std::vector< TH1 * > m_scale
unsigned long m_eventCount

Detailed Description

Definition at line 24 of file ReconDataHistogram.h.


Member Enumeration Documentation

enum ReconDataHistogram::Fig_t [private]
Enumerator:
kEnergyVsR 
kEnergyVsZ 
kRateVsR 
kRateVsZ 
kEnergy 
kZVsR 
kYVsX 

Definition at line 51 of file ReconDataHistogram.h.


Constructor & Destructor Documentation

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

Definition at line 31 of file ReconDataHistogram.cc.

  : GaudiAlgorithm(name,pSvcLocator),
    m_statsSvc(0),
    m_eventCount(0),
    m_firstTriggerTime(0)
{
}
ReconDataHistogram::~ReconDataHistogram ( ) [virtual]

Definition at line 40 of file ReconDataHistogram.cc.

{
}

Member Function Documentation

StatusCode ReconDataHistogram::initialize ( ) [virtual]

Definition at line 44 of file ReconDataHistogram.cc.

{
  // Initialize the necessary services
  StatusCode sc = service("StatisticsSvc",m_statsSvc,true);
  if(sc.isFailure()){
    error() << "Failed to get StatisticsSvc" << endreq;
    return sc;
  }
  return sc;
}
StatusCode ReconDataHistogram::finalize ( ) [virtual]

Definition at line 121 of file ReconDataHistogram.cc.

{
  info() << "Finalizing " << endreq;

  for(unsigned int Idx=0; Idx<m_scale.size(); Idx++){
    TH1* hist = m_scale[Idx];
    //hist->Sumw2();
    double duration = m_lastTriggerTime.GetSeconds() - m_firstTriggerTime->GetSeconds();
    //hist->Scale(1/duration);
    hist->SetBinContent(0,duration);
  }
  m_scale.clear();

  if( m_statsSvc ) {
    m_statsSvc->release();
  }

  if(m_firstTriggerTime){
    delete m_firstTriggerTime;
    m_firstTriggerTime = 0;
  }

  return StatusCode::SUCCESS;
}
StatusCode ReconDataHistogram::execute ( ) [virtual]

Definition at line 55 of file ReconDataHistogram.cc.

{
  verbose() << "executing: " << m_eventCount << endreq;
  // Add the current event into histograms
  DayaBay::RecHeader* recHeader = 
    get<DayaBay::RecHeader>("/Event/Rec/AdSimple");
  if(!recHeader){
    error() << "Failed to get reconstructed header:" << "/Event/Rec/AdSimple"
            << endreq;
    return StatusCode::FAILURE;
  }

  const DayaBay::RecTrigger* recTrigger = &(recHeader->recTrigger());
  if(!recTrigger){
    error() << "Failed to get reconstructed trigger: " << "Event/Rec/AdSimple"
            << endreq;
    return StatusCode::FAILURE;
  }
  
  if(!recTrigger->detector().isAD()){
    // Not an AD, continue
    return StatusCode::SUCCESS;
  }
 
  const DayaBay::Detector& detector = recTrigger->detector();
  info() << "Detector = " << detector << endreq;

  if(!m_firstTriggerTime){
    m_firstTriggerTime = new TimeStamp(recTrigger->triggerTime());
  }
  
  TH2F* h_energyVsR = dynamic_cast<TH2F*>(getOrMakeHist(1, detector, kEnergyVsR));
  TH2F* h_energyVsZ = dynamic_cast<TH2F*>(getOrMakeHist(1, detector, kEnergyVsZ));
  TH1F* h_rateVsR = dynamic_cast<TH1F*>(getOrMakeHist(1, detector, kRateVsR));
  TH1F* h_rateVsZ = dynamic_cast<TH1F*>(getOrMakeHist(1, detector, kRateVsZ));
  TH1F* h_energy = dynamic_cast<TH1F*>(getOrMakeHist(1, detector, kEnergy));
  TH2F* h_ZVsR = dynamic_cast<TH2F*>(getOrMakeHist(1, detector, kZVsR));
  TH2F* h_YVsX = dynamic_cast<TH2F*>(getOrMakeHist(1, detector, kYVsX));
  
  if(recTrigger->positionStatus()==ReconStatus::kGood){
    double x = recTrigger->position().x() / 1000;
    double y = recTrigger->position().y() / 1000;
    double z = recTrigger->position().z() / 1000;
    double r = TMath::Sqrt(x*x+y*y);
    h_rateVsR->Fill(r);
    h_rateVsZ->Fill(z);
    h_ZVsR->Fill(r,z);
    h_YVsX->Fill(x,y);
    if(recTrigger->energyStatus()==ReconStatus::kGood){
      double energy = recTrigger->energy();
      h_energyVsR->Fill(r, energy);
      h_energyVsZ->Fill(z,energy);
      h_energy->Fill(energy);
    }
  }

  m_lastTriggerTime = recTrigger->triggerTime();
  
  ++m_eventCount;
  if (0 == (m_eventCount % 1000)) {
    info() << " Processed " << m_eventCount << " events" << endreq;
  }

  return StatusCode::SUCCESS;
}
TH1 * ReconDataHistogram::getOrMakeHist ( int  run,
const DayaBay::Detector detector,
int  histIndex 
) [private]

Definition at line 405 of file DaqDataHistogram.cc.

{
  debug() << "Detector=" << detector << "," << detector.fullPackedData() << endreq;

  int detectorInt = detector.fullPackedData();
  std::map<int,TH1**>::iterator histIter = m_hist.find(detectorInt);
  TH1** histograms = 0;
  if( histIter == m_hist.end() ){
    info() << "Create " << MAXHISTS << " histograms for " << detector << endreq;
    // Initialize histogram 
    histograms = new TH1*[MAXHISTS];
    for(int i=0; i<MAXHISTS; i++) histograms[i] = 0;
    m_hist[detectorInt] = histograms;
  }else{
    // Found detector
    histograms = histIter->second;
  }

  int iHist = board * NCONNECTORS * NHISTOGRAMS + connector * NHISTOGRAMS + histIndex;
  debug() << "Histogram index = " << iHist << endreq;
  TH1* hist = histograms[board * NCONNECTORS * NHISTOGRAMS + connector * NHISTOGRAMS + histIndex];
  if(!hist){
    // Make this histogram
    std::string histName;
    if(detector.detectorId()){
      if(board){
        // Make channel histograms
        switch(histIndex){
          case kTdc:
            hist = new TH1I("TDC","TDC",1000,300,1300);
            hist->GetXaxis()->SetTitle("TDC");
            hist->GetYaxis()->SetTitle("Entries");
            break;
          case kAdcFine:
            // ADC Values (Fine Range)
            hist = new TH1I("ADCFine","ADC (Fine Range)",4096,0,4096);
            hist->GetXaxis()->SetTitle("ADC Fine");
            hist->GetYaxis()->SetTitle("Entries");
            break;
          case kAdcCoarse:
            // ADC Values (Coarse Range)
            hist = new TH1I("ADCCoarse","ADC (Coarse Range)",4096,0,4096);
            hist->GetXaxis()->SetTitle("ADC Coarse");
            hist->GetYaxis()->SetTitle("Entries");
            break;
          case kPeakCycle:
            hist = new TH1I("PeakCycle","PeakCycle",14,0,14);
            hist->GetXaxis()->SetTitle("PeakCycle");
            hist->GetYaxis()->SetTitle("Entries");
            break;
          case kRange:
            hist = new TH1I("Range","Range",3,0,3);
            hist->GetXaxis()->SetTitle("Range");
            hist->GetYaxis()->SetTitle("Entries");
            break;
          case kHitCount:
            hist = new TH1I("HitCount","HitCount",20,0,20);
            hist->GetXaxis()->SetTitle("HitCount");
            hist->GetYaxis()->SetTitle("Entries");
            break;
          case kPreAdcFine:
            // PreADC Values (Fine Range)
            hist = new TH1I("PreADCFine","PreADC (Fine Range)",200,0,200);
            hist->GetXaxis()->SetTitle("PreADC Fine");
            hist->GetYaxis()->SetTitle("Entries");
            break;
          case kPreAdcCoarse:
            // Pre ADC Values (Coarse Range)
            hist = new TH1I("PreADCCoarse","PreADC (Coarse Range)",200,0,200);
            hist->GetXaxis()->SetTitle("PreADC Coarse");
            hist->GetYaxis()->SetTitle("Entries");
            break;
          case kDarkNoise:
            hist = new TH1I("DarkNoise","DarkNoise",250,-50,200);
            hist->GetXaxis()->SetTitle("Dark Noise");
            hist->GetYaxis()->SetTitle("Entries");
            break;



        case kRandTrig:
            hist = new TH1I("RandTrig","RandTrig",250,-50,200);
            hist->GetXaxis()->SetTitle("Random Trigger");
            hist->GetYaxis()->SetTitle("Entries");
            break;



          //case kDarkNoiseVsTimePMT:
          //  hist = new TH1F("DarkNoiseVsTimePMT", "DarkNoiseVsTimePMT",
          //      60/30, m_firstTriggerTime->GetSec(), m_firstTriggerTime->GetSec() + 480);
          //  hist->GetXaxis()->SetTitle("Run Time");
          //  hist->GetYaxis()->SetTitle("Dark Noise [kHz]");
          //  hist->GetXaxis()->SetLabelSize(0.03);
          //  hist->GetXaxis()->SetTimeDisplay(1);
          //  hist->GetXaxis()->SetTimeFormat(timeFormat);
          //  hist->GetXaxis()->SetTimeOffset(0,"gmt");
          //  break;
          default:
            error() << "Unknown histogram for single channel: " << histIndex << endreq;
            return 0;
        }
      } else {
        // Make detector histIndex
        switch(histIndex){
          case kTriggerType:
            hist = new TH1F("h_sum_triggerType", "h_triggerType;triggerType", 16, 0, 16);
            for(int i = 0; i < 16; i++){
              hist->GetXaxis()->SetBinLabel(i + 1, m_pmtTrigType[i]);
            }
            break;
          //case kRateVsTriggerType:
          //  hist = new TH1F("h_sum_rateVsTriggerType", "h_rateVsTriggerType;triggerType", 16, 0, 16);
          //  hist->GetYaxis()->SetTitle("Rate");
          //  for(int i = 0; i < 16; i++){
          //    hist->GetXaxis()->SetBinLabel(i + 1, m_pmtTrigType[i]);
          //  }
          //  m_scale = hist; 
          //  //m_scale.push_back(hist); 
          //  break;
          case kAdcSum:
            hist = new TH1F("h_sum_FEEADCSUM", "h_FEEADCSUM", 1000 , 0, m_adcSumMax);
            break;
          case kNchannel:
            hist = new TH1F("h_sum_nChannel","Number of Channels in Event",200,0,200);
            hist->GetXaxis()->SetTitle("Number of Channels");
            break;
          case kHitRate:
            hist = new TH2F("h_sum_hitRate", "hitRate", 
                20,0.5,20.5,16,0.5,16.5);
            hist->GetXaxis()->SetTitle("FEE Board");
            hist->GetYaxis()->SetTitle("FEE Connector");
            break;
          case kHitNumber:
            hist = new TH1F("h_tmp_hitNumber", "h_hitNumber;channelIndex;hit_number", 
                MAXMODNUM*MAXCHNNUM , 0.5, MAXMODNUM*MAXCHNNUM+0.5);
            break;
          case kAdcMean:
            hist = new TH2F("h_sum_ADCMean", "Mean Fine ADC", 
                20,0.5,20.5,16,0.5,16.5);
            hist->GetXaxis()->SetTitle("FEE Board");
            hist->GetYaxis()->SetTitle("FEE Connector");
            break;
          case kAdcRMS:
            hist = new TH2F("h_sum_ADCRMS", "Fine ADC RMS", 
                20,0.5,20.5,16,0.5,16.5);
            hist->GetXaxis()->SetTitle("FEE Board");
            hist->GetYaxis()->SetTitle("FEE Connector");
            break;
          case kTdcMean:
            hist = new TH2F("h_sum_TDCMean", "TDC Mean", 
                20,0.5,20.5,16,0.5,16.5);
            hist->GetXaxis()->SetTitle("FEE Board");
            hist->GetYaxis()->SetTitle("FEE Connector");
            break;
          case kTdcRMS:
            hist = new TH2F("h_sum_TDCRMS", "TDC RMS", 
                20,0.5,20.5,16,0.5,16.5);
            hist->GetXaxis()->SetTitle("FEE Board");
            hist->GetYaxis()->SetTitle("FEE Connector");
            break;
          case kPreAdcMean:
            hist = new TH2F("h_sum_preADCMean", "h_preADCMean", 
                20,0.5,20.5,16,0.5,16.5);
            hist->GetXaxis()->SetTitle("FEE Board");
            hist->GetYaxis()->SetTitle("FEE Connector");
            break;
          case kPreAdcRMS:
            hist = new TH2F("h_sum_preADCRMS", "preADC RMS", 
                20,0.5,20.5,16,0.5,16.5);
            hist->GetXaxis()->SetTitle("FEE Board");
            hist->GetYaxis()->SetTitle("FEE Connector");
            break;
          case kAdcVsChannel:
            hist = new TH2F("h_sum_ADC_Channel", "h_sum_ADC_Channel", 
                MAXMODNUM*MAXCHNNUM , 0.5, MAXMODNUM*MAXCHNNUM+0.5, 
                (int)(4096*m_highGainFactor), 0, 4096*m_highGainFactor);
            hist->GetXaxis()->SetTitle("Board#times16+Channel");
            hist->GetYaxis()->SetTitle("ADC");
            break;
          case kTdcVsChannel:
            hist = new TH2F("h_sum_TDC_Channel", "h_sum_TDC_Channel", 
                MAXMODNUM*MAXCHNNUM , 0.5, MAXMODNUM*MAXCHNNUM+0.5, 1000, 300, 1300);
            hist->GetXaxis()->SetTitle("Board#times16+Channel");
            hist->GetYaxis()->SetTitle("TDC");
            break;
          case kPreAdcVsChannel:
            hist = new TH2F("h_sum_preADC_Channel", "h_sum_preADC_Channel", 
                MAXMODNUM*MAXCHNNUM , 0.5, MAXMODNUM*MAXCHNNUM+0.5, 
                (int)(4096*m_highGainFactor), 0, 4096*m_highGainFactor);
            hist->GetXaxis()->SetTitle("Board#times16+Channel");
            hist->GetYaxis()->SetTitle("preADC");
            break;
          case kDarkRate:
            hist = new TH1F("h_sum_DarkRate", "h_DarkRate", 
                MAXMODNUM*MAXCHNNUM , 0.5, MAXMODNUM*MAXCHNNUM+0.5);
            hist->GetXaxis()->SetTitle("Board#times16+Channel");
            hist->GetYaxis()->SetTitle("Dark Rate (kHz)");
            break;
          case kDarkNoiseVsChannel:
            hist = new TProfile("h_sum_DarkNoise","h_sum_DarkNoise",
                MAXMODNUM*MAXCHNNUM , 0.5, MAXMODNUM*MAXCHNNUM+0.5,-50,250);
            hist->GetXaxis()->SetTitle("Board#times16+Channel");            
            hist->GetYaxis()->SetTitle("Dark Noise (ADC)");
            break;
          case kAdcSumVsTriggerType:
            hist = new TH2F("h_sum_ADCSum_TriggerType","h_sum_ADCSum_TriggerType", 
                16, 0, 16, 1000 , 0, m_adcSumMax);
            hist->GetXaxis()->SetTitle("Trigger Type");
            for(int i = 0; i < 16; i++){
              hist->GetXaxis()->SetBinLabel(i + 1, m_pmtTrigType[i]);
            }
            hist->GetYaxis()->SetTitle("ADCSum");
            break;
          case kTriggerBlocked:
            hist = new TH1F("h_sum_TriggerBlocked", "h_sum_TriggerBlocked",
                60, m_firstTriggerTime->GetSec(), m_firstTriggerTime->GetSec() + 60 ) ;
            hist->GetXaxis()->SetTitle("Run Time (Beijing)");
            hist->GetYaxis()->SetTitle("nBlocked Triggers");
            hist->GetXaxis()->SetLabelSize(0.03);
            hist->GetXaxis()->SetTimeDisplay(1);
            hist->GetXaxis()->SetTimeFormat(timeFormat);
            hist->GetXaxis()->SetTimeOffset(0,"gmt");
            hist->GetXaxis()->SetNdivisions(505);
            break;
          //case kGain:
          //  hist = new TH1F("h_sum_Gain", "Gain vs. Time (Based on random trigger)",
          //      60, m_firstTriggerTime->GetSec(), m_firstTriggerTime->GetSec() + 60);
          //  hist->GetXaxis()->SetTitle("Run time");
          //  hist->GetYaxis()->SetTitle("PMT gain");
          //  hist->GetXaxis()->SetLabelSize(0.03);
          //  hist->GetXaxis()->SetTimeDisplay(1);
          //  hist->GetXaxis()->SetTimeFormat(timeFormat);
          //  hist->GetXaxis()->SetTimeOffset(0,"gmt");
          //  hist->GetXaxis()->SetNdivisions(505);
          //  break;
          //case kDarkNoiseVsTime:
          //  hist = new TH1F("h_sum_DarkNoiseVsTime", "DarkNoise vs. Time (Based on random trigger)",
          //      60/30, m_firstTriggerTime->GetSec(), m_firstTriggerTime->GetSec() + 480);
          //  hist->GetXaxis()->SetTitle("Run time");
          //  hist->GetYaxis()->SetTitle("PMT dark noise [kHz]");
          //  hist->GetXaxis()->SetLabelSize(0.03);
          //  hist->GetXaxis()->SetTimeDisplay(1);
          //  hist->GetXaxis()->SetTimeFormat(timeFormat);
          //  hist->GetXaxis()->SetTimeOffset(0,"gmt");
          //  hist->GetXaxis()->SetNdivisions(505);
          //  break;
          default:
            error() << "Unknown histogram for event: " << histIndex << endreq;
            return 0;
        }
      }
    }
    debug() << "Making histogram: " << hist->GetName() << endreq;
    m_statsSvc->put( getPath(run, detector, board, connector, hist->GetName()), hist );
    histograms[ board * NCONNECTORS * NHISTOGRAMS + connector * NHISTOGRAMS + histIndex] = hist;
  }

  return hist;
}
std::string ReconDataHistogram::getPath ( int  run,
const DayaBay::Detector detector,
const char *  histName 
) [private]

Definition at line 244 of file ReconDataHistogram.cc.

{
  // Construct histogram path in statistics service
  std::stringstream path;
  path << "/file1/";
  path << Site::AsString(detector.site()) << "/" << DetectorId::AsString(detector.detectorId());
  path << "/Recon";
  path << "/" << histName;
  debug() << "Histogram path = " << path << endreq;
  return path.str();
}

Member Data Documentation

IStatisticsSvc* ReconDataHistogram::m_statsSvc [private]

Definition at line 36 of file ReconDataHistogram.h.

Definition at line 37 of file ReconDataHistogram.h.

Definition at line 38 of file ReconDataHistogram.h.

std::map<int,TH1**> ReconDataHistogram::m_hist [private]

Definition at line 44 of file ReconDataHistogram.h.

std::vector<TH1*> ReconDataHistogram::m_scale [private]

Definition at line 45 of file ReconDataHistogram.h.

unsigned long ReconDataHistogram::m_eventCount [private]

Definition at line 49 of file ReconDataHistogram.h.


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:02:12 for DQMRawData by doxygen 1.7.4