GSI/ADFWatchers.C

Go to the documentation of this file.

/***************************************************************************
 *   Copyright (C) 2010 by Olivier Stezowski                               *
 *   stezow(AT)ipnl.in2p3.fr                                               *
 *                                                                         *
 *   This program is free software; you can redistribute it and/or modify  *
 *   it under the terms of the GNU General Public License as published by  *
 *   the Free Software Foundation; either version 2 of the License, or     *
 *   (at your option) any later version.                                   *
 *                                                                         *
 *   This program is distributed in the hope that it will be useful,       *
 *   but WITHOUT ANY WARRANTY; without even the implied warranty of        *
 *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the         *
 *   GNU General Public License for more details.                          *
 *                                                                         *
 *   You should have received a copy of the GNU General Public License     *
 *   along with this program; if not, write to the                         *
 *   Free Software Foundation, Inc.,                                       *
 *   59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.             *
 ***************************************************************************/

// Watcher definition
#ifndef _ADFWatchers
#include "ADFWatchers.h"
#endif

#include "TSystem.h"
#include <stdexcept>

using namespace ADF;

ClassImp(ADFWatcher);

Bool_t ADFWatcher::SetTrigger( ADF::DFTrigger *trigger ) 
{
        if ( ! Watcher::SetTrigger(trigger) ){ 
                fFrame = 0x0; 
                return false; 
        }

        // set the unique frame to be watched
        fFrame = GetTrigger()->GetInputSharedFP();
        
        return true;    
}

ClassImp(Dumper);
void Dumper::Exec(Option_t */*option*/) 
{
// Frame by Frame action
        // be sure the key (the frame) has been set properly
        if ( fFrame == 0x0 ) 
                return;
                
        // if new file is asked, next open erase the file  
        if ( fIsNew )
                { fFrame->GetFrame()->Dump(fFilename.Data(), fIsNew); fIsNew = false; }
        else
                { fFrame->GetFrame()->Dump(fFilename.Data(), fIsNew); }
}
void Dumper::SetOutputFile(const Char_t *filename)
{
        fFilename = filename; fIsNew = true;
}

ClassImp(AgataKeyWatcher) ;
        AgataKeyWatcher::AgataKeyWatcher(const Char_t *name, const Char_t *title, TDirectory *sp_dir, TDirectory *tag_dir): 
        ADFWatcher(name,title,sp_dir,tag_dir),
        fKey(0x0),
        fWhichBin(-1),
        fKeyFilter("*",true),
        fSubKeyFilter("*",true),
        fLastEventNumber(0u),
        flastTimestamp(0UL),
    fShowGlobalConf(false),
        fLengthFrame(0x0),
        fKeyMessage(0x0),
        fEventNumber(0x0),
        fTimeStamp(0x0),
        fCEventNumber(0x0),
        fCTimeStamp(0x0),
        fCoincTimeStamp(0x0),
        fLengthSubFrame(0x0),
        fKeyMessageSubFrame(0x0)
{
// allocation of the histograms you would like to fill 
        fLengthFrame = MakeTH1<TH1F>("LengthFrame","Length of the Frames (bytes)", aKByte,0,20*aKByte);
        fKeyMessage  = MakeTH1<TH1F>("KeyMessage","Message type", 3,0,3) ;
        fKeyMessage->SetBit(TH1::kCanRebin);
        fEventNumber = MakeTH1<TH1D>("EventNumber","Last event numbers", kMaxCycle,0.,kMaxCycle) ;
        fTimeStamp   = MakeTH1<TH1D>("TimeStamp","Last time stamps", kMaxCycle,0.,kMaxCycle) ;
        fCEventNumber = MakeTH1<TH1F>("CEventNumber","Consecutive event numbers", 1000,0, 10000) ;
        fCTimeStamp = MakeTH1<TH1F>("CTimeStamp","Consecutive time stamps", 1000,0, 10000000) ;
        fCoincTimeStamp = MakeTH1<TH1F>("CoincTimeStamp","Time stamp distribution in composite frames", 500,0,1000 ) ;
        fLengthSubFrame = MakeTH1<TH1F>("LengthSubFrame","Length of the SubFrames (bytes)", aKByte,0,20*aKByte);
        fKeyMessageSubFrame = MakeTH1<TH1F>("KeyMessageSubFrame","Message type for SubFrames in Composite Frame", 1,0,1);
        fKeyMessageSubFrame->SetBit(TH1::kCanRebin);
}

AgataKeyWatcher::~AgataKeyWatcher() 
{
// destruction of the histograms done because they have been added to the pool
}

void AgataKeyWatcher::SetKeyFilter(const Char_t *keypattern, const Char_t *subkeypattern)
{
    TString stmp1(keypattern), stmp2(subkeypattern);
    
    Bool_t wildcard1 = false, wildcard2 = false;
    if ( stmp1.Contains("*") )
        wildcard1 = true;
    if ( stmp2.Contains("*") )
        wildcard2 = true;
    
        TRegexp tmp1(keypattern,wildcard1), tmp2(subkeypattern,wildcard2);
        
        // change the filters only if it is a valid expression
        if ( tmp1.Status() == TRegexp::kOK ) {
                fKeyFilter = tmp1;
    }
        if ( tmp2.Status() == TRegexp::kOK ) {
                fSubKeyFilter = tmp2;
    }
}

void AgataKeyWatcher::DoCanvas(TCanvas *c, Option_t *o)
{
        TString option = o;
        if ( option == "" ) {
                c->Divide(3,2,0.001,0.001);
                c->cd(1); 
                fLengthFrame->Draw();   
                c->cd(2); 
                fLengthSubFrame->Draw();        
                c->cd(3); 
                fEventNumber->Draw();
                c->cd(4); 
                fKeyMessage->Draw();
                c->cd(5);       
                fKeyMessageSubFrame->Draw();    
                c->cd(6); 
                fTimeStamp->Draw();                     
        }
        if ( option == "q" ) {
                c->Divide(1,3,0.001,0.001);
                c->cd(1); 
                fCEventNumber->Draw();
                c->cd(2); 
                fCTimeStamp->Draw();
                c->cd(3); 
                fCoincTimeStamp->Draw();                
        }
}
        
void AgataKeyWatcher::Exec(Option_t */*option*/) 
{
        if ( gDebug > 2 ) 
                printf("AgataKeyWatcher::Exec is called \n");
        
        UInt_t length, current, eventnumber, message; ULong64_t timestamp; TString mess, mess_sec, mess_tmp;

        if ( fFrame != 0x0 ) {
                if ( fFrame->IsValid() )
                        fKey = dynamic_cast<const AgataKey *>(fFrame->GetFrame()->GetKey());
                else 
                        fKey = 0x0;
        }
        else 
                fKey = 0x0;
        
        // be sure the key (the frame) has been set properly
        if ( fKey == 0x0 ) 
                return;
        
        // first test if a filter is applied 
        mess = AgataKeyFactory::theFactory()->GetMessage(fKey->GetMessage());
        if ( ! mess.Contains(fKeyFilter) )
                return;
                
        // now fill histograms
        fLengthFrame->Fill(fKey->GetFrameLength()); fKeyMessage->Fill(mess.Data(),1);
        if ( mess.Contains("UnknownKey") ) { // an unknown key
                std::cout << " WARNING, unknown key " << std::endl; 
                fKey->Print();
        }       
        // try and read the possible event# and timestamp (may not exists in case of low version number of agata key
        const BufferIO *buf = fFrame->GetFrame()->GetRealBuffer();
        buf->SetOffset();
        if ( fKey->GetFrameLength() < 12u ) {
                eventnumber = 0u;
                timestamp   = 0UL;
        }
        else {
                (*buf) >> eventnumber;
                if ( fKey->GetFrameLength() < 20u )
                        timestamp = 0UL;        
                else (*buf) >> timestamp;               
        }
        
        // which is negative if this method has never been called
        if ( fWhichBin == -1 ) {
                fLastEventNumber = eventnumber;
                flastTimestamp   = timestamp;
                fWhichBin = 0;
        }
        
        // fill timestamp and event number related  
        fEventNumber->Fill(fWhichBin,eventnumber);
        fTimeStamp->Fill(fWhichBin,timestamp);  
        
        fCEventNumber->Fill(eventnumber-fLastEventNumber);
        fCTimeStamp->Fill(timestamp-flastTimestamp);
        
        /*
        if ( eventnumber < fLastEventNumber ) {
                fEventNumberQuality->Fill(-1);
        //      cout << " Current event # " << eventnumber << " Previous event # " << fLastEventNumber << endl;
        }
        else 
                eventnumber == fLastEventNumber ? fEventNumberQuality->Fill(0) : fEventNumberQuality->Fill(1);
        
        if ( timestamp < flastTimestamp ) {
                fTimeStampQuality->Fill(-1);
        //      cout << " Current timestamp # " << timestamp << " Previous timestamp # " << flastTimestamp << endl;
        }
        else    
                timestamp == flastTimestamp ? fTimeStampQuality->Fill(0) : fTimeStampQuality->Fill(1); */
        
        fLastEventNumber = eventnumber;
        flastTimestamp   = timestamp;
        
        // check if message contains conf. if yes dump the content in a file
    if ( mess.Contains("conf") && fShowGlobalConf ) {
                FILE *fp; Short_t which = 0;
                while ( which < 20 ) {
                        if ( gSystem->AccessPathName(Form("./%s%d.conf",mess.Data(),which)) )
                                break;
                        which++;
                }
                fp = ::fopen(Form("./%s%d.conf",mess.Data(),which),"w"); fFrame->GetFrame()->Dump(fp); ::fclose(fp);
        }
        if ( mess.Contains("event") ) {
                
                fKeyMessageSubFrame->Fill(mess.Data(),1);
                
                // try to look for a true message type at the beginning of the buffer
                Bool_t is_a_start = false; Short_t max_loop = 0u;
                buf->SetOffset();
                while ( buf->FreeSize() > 20u ) {
                        
                        // up to ten embedded frames in one frame
                        if ( max_loop++ == 10u )
                                { fLog << " Max has been reached in scan of composite frame " << nline; break; }
                                
                        current = buf->Offset(); 
                        
                        (*buf) >> length >> message >> eventnumber >> timestamp ;
                        
                        mess_tmp = AgataKeyFactory::theFactory()->GetMessage(message);
                        
//                      cout << "[" << mess_tmp.Data() << "]" << endl;
                        if ( mess_tmp.Contains("UnknownKey") ) { // an unknown key
                                is_a_start = false;
                        }
                        else is_a_start = true;
                        
                        if ( is_a_start )  { // a known key has been found
                                
                                buf->SetOffset( current + length );
                                
                                mess_sec  = "|"; 
                                mess_sec += mess_tmp; 
                                mess_sec += "|";
                                
                                mess_sec = Form("|%s|",mess_tmp.Data());
                                
                                if ( mess_tmp.Contains(fSubKeyFilter ) ) {

                                        fLengthSubFrame->Fill(length);
                                        fKeyMessageSubFrame->Fill(mess_sec.Data(),1);   
                                        
                                        fCoincTimeStamp->Fill(timestamp-flastTimestamp);
                                }
                        }
                        else buf->SetOffset( current + 4u );
                }
        }       
        
        if ( fWhichBin == kMaxCycle ) {
                fWhichBin = 0u;
                fEventNumber->Reset("ICE"); fTimeStamp->Reset("ICE");
        }
        else fWhichBin++;
}

ClassImp(TSRate);

TSRate::TSRate(const Char_t *name, const Char_t *title, TDirectory *sp_dir, TDirectory *tag_dir) :
        ADFWatcher(name,title,sp_dir,tag_dir),
        fTSRates(0)
{ 
}

void TSRate::Exec(Option_t * /*option*/)
{       
        // fLog << debug << " In TSRate::Exec() " << nline;

        if ( fFrame && fFrame->IsValid() ) {
                AgataKey *key = dynamic_cast<AgataKey *>(fFrame->GetFrame()->GetKey());
                if ( key ) {
                        ULong64_t ts = key->GetTimeStamp() ;
                        if ( fTSRates.size() > 0 ) {
                                if ( fTSRates.back() < ts ) 
                                        fTSRates.push_back( ts ) ;
  //                              else 
 //                                     std::cout << " Two consecutive TS in a bunch inversed !" << std::endl;                                    
                        }
                        else fTSRates.push_back( ts ) ;
                }
        }
        
        // fLog << debug << " Out TSRate::Exec() " << dolog;
        
}

ClassImp(TSRateCollector);

TSRateCollector::TSRateCollector(const Char_t *name, const Char_t *title, TDirectory *sp_dir, TDirectory *tag_dir) :
        Watcher(name,title,sp_dir,tag_dir),
        fDefinitionRate(),
        fCurrentRate(),
        fTSRates(),
        fRates(),
        fPacketWidth(),
        fDeadTime(),
        fT0(kMaxULong64)
{
        fTSRates.SetOwner(false); fPacketWidth.SetOwner(false); fDeadTime.SetOwner(false);
        
        // rate every 30 second for one hour
        fDefinitionRate.Add(30*fgkOneSecond,120);
        // rate every 10 minutes for 12 hour
        fDefinitionRate.Add(600*fgkOneSecond,72);       
        // rate every hour for 48 hours
        fDefinitionRate.Add(3600*fgkOneSecond,48);
}

Bool_t TSRateCollector::SetT0()
{
        ULong64_t val, min = kMaxULong64;
                
        // try and see if first time it is called or if timestamp has been reset
        TIter next(&fTSRates);
        TSRate *ts_rate;
        while ( (ts_rate = (TSRate *)next()) ) {
                if ( ts_rate->LastValues().size() ) {
                        val = ts_rate->LastValues().at(0);      
                        if ( val < min )
                                min = val;
                }
        }
        if ( min < fT0  ) {
                // change histograms title to include starting time
                TDatime dt;
                for (Int_t i = 0; i < fRates.GetHists()->GetSize(); i++) {
                        
                        // reset stat
                        fCurrentRate[i].Reset(fDefinitionRate.fNbScales);
                        
                        if ( fRates.GetHists()->At(i) == 0x0 ) 
                                continue;
                        TString tmp = fRates.GetHists()->At(i)->GetName();
                        tmp += " since ";
                        tmp += dt.AsString();
                        
                        ((TH1 *)fRates.GetHists()->At(i))->Reset();
                        ((TH1 *)fRates.GetHists()->At(i))->SetTitle(tmp.Data());
                        
                }       
                fT0 = min;
                std::cout << " New TIMESTAMP in TSRateCollector" << fT0 << std::endl;
        }
        
        return (fT0 != kMaxULong64);
}

void TSRateCollector::AddHook(Watcher *ts_rate, const char *fd_name)
{
        const Short_t wheel_size = 5; 
        Color_t color_wheel[wheel_size] = { kBlue, kGreen+2, kOrange, kRed, kViolet }; 
        Style_t style_wheel[wheel_size] = { 3001, 3002, 3003, 3004, 3005 };

        if ( !ts_rate->InheritsFrom("TSRate") ) 
                return;
        
        if ( ts_rate ) {
                // keep an hook on it and creates the associated struct to compute stats
                fTSRates.Add ( ts_rate );
                C_Rate c_rate(fDefinitionRate.fNbScales); fCurrentRate.push_back( c_rate );
                
                // now add histograms
                TString hname = fd_name;
                hname += "_Rate";
                
                Int_t hsize = 0;
                for (Int_t i = 0; i < fDefinitionRate.fNbScales; i++)
                        hsize += fDefinitionRate.fNbSamples[i] + gEmptyBins;
                        
                hsize += 10;
                TH1F *h = MakeTH1<TH1F>(hname.Data(),"Rate",hsize,0,hsize);
                fRates.Add(h);
                h->SetFillColor(color_wheel[(fRates.GetHists()->GetSize()-1)%wheel_size]);
                h->SetFillStyle(style_wheel[(fRates.GetHists()->GetSize()-1)%wheel_size]);
                h->SetXTitle("[0][30s/bin = 1hour] -- [170][10m/bin = 12hours] -- [292][1h/bin = 2days]");
                h->SetYTitle("k/s");

                hname = fd_name; hname += "_Width";
                h = MakeTH1<TH1F>(hname.Data(),"Width distribution of packets (s)",200,0,20);
                fPacketWidth.Add(h);
                h->SetFillColor(color_wheel[(fRates.GetHists()->GetSize()-1)%wheel_size]);
                h->SetFillStyle(style_wheel[(fRates.GetHists()->GetSize()-1)%wheel_size]);
                
                hname = fd_name; hname += "_Deadtime";
                h = MakeTH1<TH1F>(hname.Data(),"Dead time distribution between consecutive packets (s)",100,0,100);
                fDeadTime.Add(h);
                h->SetFillColor(color_wheel[(fRates.GetHists()->GetSize()-1)%wheel_size]);
                h->SetFillStyle(style_wheel[(fRates.GetHists()->GetSize()-1)%wheel_size]);
        }
}

void TSRateCollector::DoCanvas(TCanvas *c, Option_t * /*o*/)
{
    Int_t nbcolumn = 3; c->Divide(nbcolumn,fRates.GetHists()->GetSize()/nbcolumn,0.001,0.001);

    for (Int_t i = 0; i < fRates.GetHists()->GetSize(); i++) {
                c->cd(i+1);
                TH1 *h = (TH1F *)fRates.GetHists()->At(i);
                if ( h )
                        h->Draw();
    }
    TString cmd = Form("RunLoupeOnPad( (TCanvas *)gROOT->GetListOfCanvases()->FindObject(\"%s\") )", c->GetName());
    gROOT->ProcessLine(cmd.Data());
}

#include "TMath.h"

void TSRateCollector::Exec(Option_t * /*option*/)
{       
        const static Double_t to_rate = fgkOneSecond/(1000.); 

        TSRate *ts_rate; TH1 *h; ULong64_t dt0, dt, slot, cycle; 

        if ( gDebug > 0 ) 
                printf("TSRateCollector::Exec is called \n");
        
        if ( fTSRates.GetSize() == 0 )
                return; // to avoid doing things in case there are no input
        
        // first time this is called, computes T0
        if ( ! SetT0() ) {
                printf(" WARNING in TSRateCollector::Exec Cannot get First timestamp \n");              
                return;
        }
        
        // loop on all registered TSRate
        
        for (Int_t i = 0; i < fTSRates.GetSize(); i++) {
                
                // get TSRate
                ts_rate = (TSRate *)fTSRates.At(i); 
                
//              std::cout << " TSRate name : " << ts_rate->GetName() 
//                              << " , registered @ position " << i << " , has view # TS : " << ts_rate->LastValues().size() << std::endl;
                
                // nothing to do
                if ( ts_rate->LastValues().size() < 2 )
                        continue;
                
                // compute some value to have an idea of the TS structure
                h = (TH1F *)fPacketWidth.At(i);
                Double_t width = ts_rate->LastValues().at(ts_rate->LastValues().size()-1) - ts_rate->LastValues().at(0);
                width /= fgkOneSecond;
                h->Fill(width);
                
                h = (TH1F *)fDeadTime.At(i);
                if ( ts_rate->LastValues().at(0) > fCurrentRate[i].fLastTS ) {
                        h->Fill((ts_rate->LastValues().at(0)-fCurrentRate[i].fLastTS)/fgkOneSecond);
                }
                else { // ste underflow
                        h->Fill(-(fCurrentRate[i].fLastTS-ts_rate->LastValues().at(0))/fgkOneSecond);           
  //                      continue;     
                } 
                // h->Fill(ts_rate->LastValues().size());
                
//              std::cout << " WIDTH " << width << std::endl;
                
                // keep last TS of the packet for next packet
                fCurrentRate[i].fLastTS = ts_rate->LastValues().at(ts_rate->LastValues().size()-1);
                
/*              // only one frame ... very low counting rate ...
                if ( ts_rate->LastValues().size() < 2 ) {
//                      ts_rate->LastValues().resize(0);
                        // do not resize so that the value is kept in the TSRate and a new one is going to be added
                        // not really a good way to proceed but should give a quite good idea in case of very low counting rates
                        continue;
                } */
                
                // statistics
                h = (TH1F *)fRates.GetHists()->At(i);
                
                for (size_t j = 1; j < ts_rate->LastValues().size(); j++) {
                        
                        if ( ts_rate->LastValues().at(j) < fT0 ) {
                                std::cout << " ERROR Timestamp in array smaller that T0 " << std::endl;
                                continue;
                        }
                                                
                        // compute time since T0 and since first timestamp
                        dt0 = ts_rate->LastValues().at(j) - fT0;
                        // compute time between two consecutive TS and the associated rate
                        dt  = ts_rate->LastValues().at(j) - ts_rate->LastValues().at(j-1); 

                        // loop on the different scales to compute statistics
                        Int_t h0 = 0;
                        for (Int_t k = 0; k < fDefinitionRate.fNbScales; k++) {
                                // cycle since T0
                                cycle = dt0 / fDefinitionRate.fPrecision[k]; 
                                //
                                if ( cycle == fCurrentRate[i].fCycle[k] ) { // add to current for integration
                                        fCurrentRate[i].fSum[k] += dt;
                                        fCurrentRate[i].fNbCounts[k]++;
                                }
                                else { // a new input in the histogram
                                        
                                        // bin number to be changed
                                        slot = fCurrentRate[i].fCycle[k] % fDefinitionRate.fNbSamples[k];
                                        
                                        // fill h
                                        // it has to be with Sum to take into account the real time counting time and remove empty cycle
/*                                      if ( slot == 0 ) {
                                                for (Int_t bin = h0; bin < h0+fDefinitionRate.fNbSamples[k]; bin++) {
                                                        h->SetBinContent(bin+1,0);
                                                }
                                        } */
                                        if ( fCurrentRate[i].fSum[k] > 0 
                                                 && ((slot+h0+1) < UInt_t(h->GetNbinsX())) 
                                                 && ((slot+h0+1) > 0) ) {      
                                                h->SetBinContent(slot+h0+1,fCurrentRate[i].fNbCounts[k]*to_rate/fCurrentRate[i].fSum[k]);
                                        }
                                        else {
                                                std::cout << " SUM is O !! " 
                                                << fCurrentRate[i].fSum[k] << " j " << j << " " << ts_rate->LastValues().size() 
                                                << " " << fCurrentRate[i].fNbCounts[k] << " " << dt << std::endl;
                                        }

/*                                      else 
                                         std::cout << slot+h0+1 << std::endl;                                   
                                         std::cout << fCurrentRate[i].fNbCounts[k]*to_rate/fCurrentRate[i].fSum[k] << std::endl; */
                                        
                                        // keep the current situation
                                        fCurrentRate[i].fSum[k] = dt;
                                        fCurrentRate[i].fNbCounts[k]  = 1;
                                        
                                        fCurrentRate[i].fCycle[k] = cycle;                                      
                                }
                                h0 += (fDefinitionRate.fNbSamples[k] + gEmptyBins);
                        }
                }
                
                ts_rate->LastValues().resize(0);
        } // i
        
        if ( 0 ) 
                printf("TSRateCollector::Exec is done \n");
}


ClassImp(TimeStampWatcher);

TimeStampWatcher::TimeStampWatcher( const char * name, const char * title, TDirectory *sp_dir, TDirectory *tag_dir):
    ADFWatcher(name,title,sp_dir,tag_dir)
{
    int TMax = 3600;
    TimeStampHist = MakeTH1<TH1F>("TimeStampHist", "Time stamp distribution",TMax*10., 0, TMax);
    FirstTimeStamp = 0;
}

TimeStampWatcher::~TimeStampWatcher()
{
}


void TimeStampWatcher::DoCanvas(TCanvas * /*c*/, Option_t *o)
{
    TimeStampHist->Draw(o);
}

void TimeStampWatcher::Exec(Option_t * /*option*/)
{

    ULong64_t fTimestamp = ((AgataKey *)fFrame->GetFrame()->GetKey())->GetTimeStamp();

    if(fTimestamp != 0 && FirstTimeStamp==0) FirstTimeStamp = fTimestamp;

    double TSinSec = (fTimestamp-FirstTimeStamp)*10*1e-9;

    TimeStampHist->Fill(TSinSec);
    LastTimeStamp = fTimestamp;

}

void TimeStampWatcher::Zero(Option_t * /*opt1*/, Option_t * /*opt2*/)
{
    TimeStampHist->Reset();
    FirstTimeStamp = LastTimeStamp;
}