/data3/calcul/jacquem/working_dir/Micromegas/micromegasFrameWork/src/analyse/root/monitoring.cpp

Go to the documentation of this file.

/* @version $Revision: 1566 $ * @modifiedby $Author: jacquem $ * @lastmodified $Date: 2012-03-08 14:00:10 +0100 (Thu, 08 Mar 2012) $ */
/* 
 * File:   beamProfile.cpp
 * Author: jacquem
 *
 * Created on May 6, 2010, 10:03 AM
 */

// System includes
#include <iostream>
#include <stdlib.h>
#include <map>
#include <string>

// Root includes
#include "TROOT.h"
#include "TGraph.h"
#include "TStyle.h"
#include <TFile.h>
#include <TTree.h>
#include <TLegend.h>
#include <TKey.h>
#include <TH1I.h>
#include <TH2I.h>
#include <TText.h>
#include <TCanvas.h>
#include <TPaveLabel.h>

// Framework inlcudes
#include "Run.hh"
#include "Detector.hh"
#include "XMLTool.hh"
#include "Toolbox.hh"
#include "MicroException.hh"
#include "Log.hh"

#include "root/MTRun.hh"
#include "root/MTEvent.hh"
#include "root/MTChannel.hh"
#include "tools/Arguments.hh"

#include "Chamber.hh"
#include "Detector.hh"
using namespace std;



int main(int argc, char** argv) {


  FILELog::ReportingLevel() = FILELog::FromString(INFO);

  Arguments arg(argc,argv,"monitoring inputFile outputFile -l evtNum");
  if (arg.getNbOfArgs() != 2 ) { arg.usage() ; exit(-1);}
  string rootName = arg.getArgument(1);
  string outputFileName = arg.getArgument(2);

  TFile *output= new TFile(outputFileName.c_str(),"RECREATE");


  Char_t flag[20];
  TTree *flagTree = new TTree("FLAG", "evant flag");
  flagTree->Branch("Flag", flag,"flag/C");
   
  vector<TH1I* > dtVec;

 map<unsigned int, map<UInt_t, Float_t>* > tempAsu1Map; // 1rst map, key is DifId , 2nd if Event Number
 map<unsigned int, map<UInt_t, Float_t>* > tempAsu2Map; // 1rst map, key is DifId , 2nd if Event Number
 map<unsigned int, map<UInt_t, Float_t>* > tempDifMap; // 1rst map, key is DifId , 2nd if Event Number


  UInt_t eventLimit = 0;
  if ( arg.getOption("-l").size() != 0 )
  {
    eventLimit = atoi(arg.getOption("-l").c_str());
  }

  // XY histos
  std::map<ui16,TH2I*> xyKeyChamber; // key is chamberId
  bool construct = false;

  TFile f(rootName.c_str());

  UInt_t nBoarderEvt = 0;
  UInt_t nOnlyFFFF = 0;
  UInt_t nOnlyAAAA = 0;
  UInt_t nOnly5555 = 0;
  UInt_t nCRCError = 0;
  UInt_t nBadChipId = 0;
  UInt_t nBadBoardId = 0;
  UInt_t nBadDifId = 0;
  UInt_t nUnknowError = 0;
  UInt_t nNotValidEvent = 0;
  TH1I validityFlag("Event_validity_flag" ,  "Event validity flag" , 5000,0.5,5000.5);

  
  
 

  TH2I* xy = NULL;
  TIter nextkey(f.GetListOfKeys());
  TKey *key;
  while (key = (TKey*)nextkey()) 
  {
    TTree *tree = (TTree*)key->ReadObj();
    FILE_LOG(logINFO  ) << "New TTree [" << tree->GetName() << endl;
    MTRun *run=NULL;
    Detector detector;
    try
    {
      UInt_t runId = 0;
      run=dynamic_cast<MTRun *>(tree->GetUserInfo()->First());
      if (construct == false)
      {
        runId = run->GetRunId();
        FILE_LOG(logINFO) << "The runId[" << runId << "]" << endl;
        if ( runId != 0)
        {
          if (! Toolbox::getDetector(runId,detector)) { exit(-1);}
          // loop over chamber config
          ChamberMap_t chambers = detector.getChambers();
          for( ChamberMap_t::const_iterator ii=chambers.begin(); ii!=chambers.end(); ++ii)
          {
            const Chamber& chamber= *((*ii).second);
            
            string xyname("hxy_");
            Toolbox::addIntToString(xyname,chamber.getId());

            string  xytitle("XY distribution of hits in chamber ");
            Toolbox::addIntToString(xytitle,chamber.getId());

            FILE_LOG(logINFO) << "Chamber id[" << chamber.getId() << "] , type: " << chamber.getType() << endl;
            // create histogram
            if ( chamber.getType().compare(MICROROCCHAMBER6) == 0 ) // MICROMEGAS
            {                                                   // 1 m -> 1 000 000 um
              xy = new TH2I(xyname.c_str(),xytitle.c_str(),193,0,965000,194,0,970000);
            }
            else if  ( chamber.getType().compare(HARDROC2IPNLCHAMBER6) == 0 ) // RPC
            {
              Int_t max = 10400*96; // 10400 micro * 96 pad
              xy = new TH2I(xyname.c_str(),xytitle.c_str(),96,0,max,96,0,max);
            }

            if ( xy != NULL )
            {
              xy->SetXTitle("y");
              xy->SetYTitle("x");
            //cout << "Create XY dist for " << xytitle << endl;

            // store it to the map
              xyKeyChamber.insert(make_pair(chamber.getId(),xy));
            }
            else { 
            FILE_LOG(logINFO) << "ELSE Chamber id[" << chamber.getId() << "] , type: " << chamber.getType() << endl;
            }
          }
        }
      construct = true;
      }
    }
    catch (MicroException e) {}

    MTEvent *evt =  new MTEvent();
    TBranch *branch= tree->GetBranch("MTEvent");
    branch->SetAddress(&evt);
    MTChannel* channel =NULL;
    int nbEntries = tree->GetEntries();
    if ( eventLimit == 0 ) { eventLimit = nbEntries; }
    unsigned int interval = (unsigned int) eventLimit / 10;
    if (interval == 0) interval = 1;
    for ( int evtNum = 0; evtNum <eventLimit ; evtNum++)
    {
      tree->GetEntry(evtNum);
 
      UInt_t nbPads=100;     // input   2/3 du nombre de pad format le cadre d'un asu 1 ASU -> ( 6 chips * 8 pads * 2 cotes + 4 chips * 8 pad * 2 cote ) / 3 * 2
      UInt_t asu;            // output
      UInt_t abs_time;       // output
      UInt_t chamber;       // output

      if ( evt->IsSquare(nbPads,asu,abs_time, chamber) == true ) { nBoarderEvt++; } 
      if ( ! evt->IsValid() ) { nNotValidEvent++; }
      if ( evt->GetCrcIsCorrect() == false ) 
      {
        nCRCError++; 
      } 
/*
#define VALID               0x1
#define UNKNOW_NOT_VALID    0x2
#define BAD_CHIP_ID         0x4
#define BAD_BOARD_ID        0x8
#define BAD_DIF_ID          0x10
#define CRC_ERROR           0x20
#define ONLY_FFFF           0x40
#define ONLY_AAAA           0x80
#define ONLY_5555           0x100
#define GLOBAL_HEADER_ERROR 0x200
#define FRAME_HEADER_ERROR  0x400
#define TRAILER_C3_FOUND    0x800
#define END_OF_FILE_ERROR   0x1000
*/
      if ( (evt->GetValidFlag() & VALID) ==  VALID )                              
      {
        validityFlag.Fill(VALID);
        strcpy(flag, "VALID") ;
        flagTree->Fill();
      }
      if ( (evt->GetValidFlag() & UNKNOW_NOT_VALID) ==  UNKNOW_NOT_VALID )        
      { 
        nUnknowError++;
        strcpy(flag, "UNKNOW_NOT_VALID") ;
        flagTree->Fill();
      }
      if ( (evt->GetValidFlag() & ONLY_AAAA) == ONLY_AAAA )                       
      {
        nOnlyAAAA++;
        validityFlag.Fill(ONLY_AAAA);
        strcpy(flag, "ONLY_AAAA") ;
        flagTree->Fill();
      } 
      if ( (evt->GetValidFlag() & ONLY_5555) == ONLY_5555 )                       
      { 
        nOnly5555++; 
        validityFlag.Fill(ONLY_5555);
        strcpy(flag, "ONLY_5555") ;
        flagTree->Fill();
      } 
      if ( (evt->GetValidFlag() & ONLY_FFFF) == ONLY_FFFF )                       
      {
         nOnlyFFFF++;validityFlag.Fill(ONLY_FFFF); 
        strcpy(flag, "ONLY_FFFF") ;
        flagTree->Fill();
      }
      if ( (evt->GetValidFlag() & BAD_CHIP_ID) ==  BAD_CHIP_ID )                  
      { 
        nBadChipId++;
        strcpy(flag, "BAD_CHIP_ID") ;
        flagTree->Fill();
        validityFlag.Fill(BAD_CHIP_ID); 
      }
      if ( (evt->GetValidFlag() & BAD_BOARD_ID) ==  BAD_BOARD_ID )                
      {
        nBadBoardId++;
        strcpy(flag, "BAD_BOARD_ID") ;
        flagTree->Fill();
        validityFlag.Fill(BAD_BOARD_ID); 
      }
      if ( (evt->GetValidFlag() & BAD_DIF_ID) ==  BAD_DIF_ID )                    
      { 
        nBadDifId++;
        strcpy(flag, "BAD_DIF_ID") ;
        flagTree->Fill();
        validityFlag.Fill(BAD_DIF_ID); 
      }
      if ( (evt->GetValidFlag() & GLOBAL_HEADER_ERROR) ==  GLOBAL_HEADER_ERROR )  
      {
        validityFlag.Fill(GLOBAL_HEADER_ERROR); 
        strcpy(flag, "GLOBAL_HEADER_ERROR") ;
        flagTree->Fill();
      }
      if ( (evt->GetValidFlag() & FRAME_HEADER_ERROR) ==  FRAME_HEADER_ERROR )    
      {
        validityFlag.Fill(FRAME_HEADER_ERROR);
        strcpy(flag, "FRAME_HEADER_ERROR") ;
        flagTree->Fill();
      }
      if ( (evt->GetValidFlag() & TRAILER_C3_FOUND) ==  TRAILER_C3_FOUND )        
      {
        validityFlag.Fill(TRAILER_C3_FOUND); 
        strcpy(flag, "TRAILER_C3_FOUND") ;
        flagTree->Fill();
      }
      if ( (evt->GetValidFlag() & END_OF_FILE_ERROR) ==  END_OF_FILE_ERROR )      
      {
        validityFlag.Fill(END_OF_FILE_ERROR); 
        strcpy(flag, "END_OF_FILE_ERROR") ;
        flagTree->Fill();
      }
      if ( (evt->GetValidFlag() & CRC_ERROR) ==  CRC_ERROR )                        
      {
        validityFlag.Fill(CRC_ERROR); 
        strcpy(flag, "CRC_ERROR") ;
        flagTree->Fill();
      }

      



      map<UInt_t, UInt_t> dtMap; // key is BcIdDif() - BcIdHit. Value is number of entry
      UInt_t maxDt = 0;
      UInt_t minDt = 0xFFFFFFFF;

      int nbChannel = evt->GetNchannel();
      for(int i=0;i<nbChannel  ;i++)
      {
        channel = (MTChannel*)evt->GetChannels()->UncheckedAt(i);
        UInt_t t2 = channel->GetBcIdDif();
        UInt_t t3 = channel->GetBcIdHit();
        UInt_t dt = t2 - t3;
        if ( dt > maxDt ) { maxDt = dt ; }
        if ( dt < minDt ) { minDt = dt ; }


        if (dtMap.find(dt) == dtMap.end())
        { //if this dt doesn't exist yet
          dtMap.insert(make_pair(dt,1));
        }
        else
        { //if this abs time already exist
          dtMap[dt] = dtMap[dt]++;
        }

  

        TH2I *hist = xyKeyChamber[channel->GetChamberId()];
        if ( hist ) 
        {
          if ( detector.getChamberById(channel->GetChamberId()).getType().compare(MICROROCCHAMBER6) == 0 ) // MICROMEGAS
          {
              hist->Fill(channel->GetY() + 2500,channel->GetX()+2500,channel->GetDigitalValue());
              hist->Fill(channel->GetY()+2500,channel->GetX()-2500,channel->GetDigitalValue());
              hist->Fill(channel->GetY()-2500,channel->GetX()+2500,channel->GetDigitalValue());
              hist->Fill(channel->GetY()-2500,channel->GetX()-2500,channel->GetDigitalValue());
          }
          else if ( detector.getChamberById(channel->GetChamberId()).getType().compare(HARDROC2IPNLCHAMBER6) == 0 ) 
          { 
              hist->Fill(channel->GetY(),channel->GetX(),channel->GetDigitalValue());
          }

        }
      }
      
      if ( evtNum % interval == 0 )
      {
        //cout << "Min Dt[" << minDt <<  "], maxDt[" << maxDt << "]" << flush;
        string name = "dt_" ;
        name = name + tree->GetName();
        Toolbox::addIntToString(name,evtNum);
        string title = "BcIdDif - BcIdHit";
        TH1I *dtHist = new TH1I(name.c_str(),title.c_str(), maxDt-minDt,  0 ,maxDt);
        
        name = "dt_projection";
        name = name + tree->GetName();
        Toolbox::addIntToString(name,evtNum);

        TH1I *dtProjHist = new TH1I(name.c_str(),title.c_str(), nbChannel,  0 ,nbChannel);
        for ( map<UInt_t ,UInt_t>::const_iterator iterC = dtMap.begin(); iterC != dtMap.end(); iterC++)
        {
          Int_t val = (*iterC).second;
          dtHist->Fill((*iterC).first,val);
          dtProjHist->Fill(val);
        }
        dtVec.push_back(dtHist);
        dtVec.push_back(dtProjHist);
      }

      // Asu1, 2 and Dif temperature
      const map<unsigned int , float>& tempAsu1EvtMap = evt->GetTemperatureAsu1();
      const map<unsigned int , float>& tempAsu2EvtMap = evt->GetTemperatureAsu2();
      const map<unsigned int , float>& tempDifEvtMap = evt->GetTemperatureDif();

      for ( map<unsigned int , float>::const_iterator it = tempAsu1EvtMap.begin() ; it != tempAsu1EvtMap.end(); it++)
      {
        unsigned int difId = it->first ;
        float asu1temp = it->second;
        float asu2temp = tempAsu2EvtMap.find(difId)->second;
        float diftemp = tempDifEvtMap.find(difId)->second;
      //if ( difId == 89 ) cout << "eventId[" << evtNum << "], asu1temp[" << asu1temp << "],asu2temp[" << asu2temp << "],difTemp[" << diftemp << "]"<< endl;

        map<UInt_t , Float_t>* foo = NULL;
        if ( tempAsu1Map.find(difId) == tempAsu1Map.end())
        {
          foo = new map<UInt_t , Float_t>();
          tempAsu1Map.insert(make_pair(difId,new map<UInt_t , Float_t>()));
          tempAsu2Map.insert(make_pair(difId,new map<UInt_t , Float_t>()));
          tempDifMap.insert(make_pair(difId,new map<UInt_t , Float_t>()));
        }
        if ( asu1temp != 0. ) (tempAsu1Map.find(difId)->second)->insert(make_pair(evtNum,asu1temp));
        if ( asu2temp != 0. ) (tempAsu2Map.find(difId)->second)->insert(make_pair(evtNum,asu2temp));
        if ( diftemp != 0. ) (tempDifMap.find(difId)->second)->insert(make_pair(evtNum,diftemp));
      }

      FILE_LOG(logINFO  ) << "Write Event[" << evtNum << "] \t\t\t\t\t\t\t\t\t\t\t\r" << flush;
    }
 }




   //TCanvas *statistic = new TCanvas("statistic", "statistic",100,10,700,500);
   TCanvas *statistic = new TCanvas("statistic", "statistic",100,10,400,300);
   statistic->Range(0,0,20,24);
   statistic->SetFillColor(10);
   statistic->SetBorderSize(4);

   TPaveLabel *pl = new TPaveLabel(3,22,17,23.7, "Monitoring result","br");
   pl->SetFillColor(18);
   pl->Draw();

   TText t(0,0,"a");
   t.SetTextFont(62);
   t.SetTextSize(0.025);
   t.SetTextAlign(12);

   string msg = "Total of not valid Event: ";
   Toolbox::addIntToString(msg,   nNotValidEvent); msg = msg + ( " / "); Toolbox::addIntToString(msg,   eventLimit );
   float p = float(float(nNotValidEvent) / float( eventLimit) * 100 );
   msg = msg + ( "-> "); Toolbox::addFloatToString(msg,  p ); msg = msg + ( "% ");
   t.DrawText(2,20,msg.c_str());

   msg = "Number of Event with CRC incorrect: ";
   Toolbox::addIntToString(msg,   nCRCError); msg = msg + ( " / "); Toolbox::addIntToString(msg,   eventLimit );
   p = float(float(nCRCError) / float( eventLimit) * 100 );
   msg = msg + ( "-> "); Toolbox::addFloatToString(msg, p ); msg = msg + ( "% ");
   t.DrawText(2,18,msg.c_str());
   
   msg = "Number of not valid Event because of  bad chip id: ";
   Toolbox::addIntToString(msg,   nBadChipId); msg = msg + ( " / "); Toolbox::addIntToString(msg,   eventLimit );
   p = float(float(nBadChipId) / float( eventLimit) * 100 );
   msg = msg + ( "-> "); Toolbox::addFloatToString(msg,  p ); msg = msg + ( "% ");
   t.DrawText(2,16,msg.c_str());

   msg = "Number of not valid Event because of unknow error: ";
   Toolbox::addIntToString(msg,   nUnknowError); msg = msg + ( " / "); Toolbox::addIntToString(msg,   eventLimit );
   p = float(float(nUnknowError) / float( eventLimit) * 100 );
   msg = msg + ( "-> "); Toolbox::addFloatToString(msg,  p ); msg = msg + ( "% ");
   t.DrawText(2,14,msg.c_str());

   msg = "Number of Event with  only FFFF as data ";
   Toolbox::addIntToString(msg,   nOnlyFFFF); msg = msg + ( " / "); Toolbox::addIntToString(msg,   eventLimit );
   p = float(float(nOnlyFFFF) / float( eventLimit) * 100 );
   msg = msg + ( "-> "); Toolbox::addFloatToString(msg,  p ); msg = msg + ( "% ");
   t.DrawText(2,12,msg.c_str());


   msg = "Number of Event with  only AAAA as data ";
   Toolbox::addIntToString(msg,   nOnlyAAAA); msg = msg + ( " / "); Toolbox::addIntToString(msg,   eventLimit );
   p = float(float(nOnlyAAAA) / float( eventLimit) * 100 );
   msg = msg + ( "-> "); Toolbox::addFloatToString(msg,  p ); msg = msg + ( "% ");
   t.DrawText(2,10,msg.c_str());

   msg = "Number of boarder Event (with  a least 100 boarder pad): ";
   Toolbox::addIntToString(msg,   nBoarderEvt); msg = msg + ( " / "); Toolbox::addIntToString(msg,   eventLimit );
   p = float(float(nBoarderEvt) / float( eventLimit) * 100 );
   msg = msg + ( "-> "); Toolbox::addFloatToString(msg,  p ); msg = msg + ( "% ");
   t.DrawText(2,8,msg.c_str());

   statistic->cd();
   output->cd();

  for ( map<ui16,TH2I*>::const_iterator iterC = xyKeyChamber.begin(); iterC != xyKeyChamber.end(); iterC++)
  {
    TH2I *hist = (*iterC).second;
    if ( hist )  
    {  
      hist->Write();
    } 
  }

 
  TCanvas deltaT("deltatT");
  double racine = sqrt(dtVec.size());
  UInt_t intSqrt = (UInt_t)(racine+1);

  deltaT.Divide(intSqrt,intSqrt);
  
  UInt_t index = 1;

  for (vector<TH1I* >::const_iterator iterC = dtVec.begin(); iterC != dtVec.end(); iterC++)
  {
    TH1I *hist = *iterC;
    if ( hist )
    { 
      if ( index % 2 == 0) deltaT.GetPad(index)->SetLogx();
      deltaT.GetPad(index)->SetLogy();
      deltaT.cd(index);
      index++;
      hist->Draw();
      hist->Write();
    }
  }

  
  for ( map<unsigned int,map< UInt_t, Float_t>*>::const_iterator iterC = tempAsu1Map.begin(); iterC != tempAsu1Map.end(); iterC++)
  {
    UInt_t difId = (*iterC).first;
    string title("temperature_dif_");
    Toolbox::addIntToString(title,difId);
    TCanvas temperature(title.c_str());
    temperature.Divide(1,3);

    map<UInt_t, Float_t>* asu1 = (*iterC).second;
    map<UInt_t, Float_t> *asu2 = tempAsu2Map.find(difId)->second;
    map<UInt_t, Float_t> *dif = tempDifMap.find(difId)->second;
    
    
    UInt_t n = asu1->size();
    Float_t xAsu1[n], yAsu1[n];
    Float_t xAsu2[n], yAsu2[n];
    Float_t xDif[n], yDif[n];
    index = 0;
    if ( asu1 )  
    {
      for ( map< UInt_t, Float_t>::const_iterator iterTemp = asu1->begin(); iterTemp != asu1->end(); iterTemp++)
      {   
        UInt_t eventId = (*iterTemp).first; 
        xAsu1[index] = eventId;
        yAsu1[index] = (*iterTemp).second;
        xAsu2[index] = eventId;
        yAsu2[index] = asu2->find(eventId)->second;
        xDif[index] = eventId;
        yDif[index] = dif->find(eventId)->second;
//      if ( difId == 89 ) cout << "yAsu1[index][" << yAsu1[index] << "],yAsu2[index][" << yAsu2[index] << "],yDif[index][" << yDif[index] << "]"<< endl;
        index++;
      }
      temperature.cd(1);
      TGraph *gr1 = new TGraph(n,xAsu1,yAsu1);
      gr1->SetTitle("Asu 1 temperature");
      gr1->Draw("AL");
      temperature.cd(2);
      TGraph *gr2 = new TGraph(n,xAsu2,yAsu2);
      gr2->SetTitle("Asu 2 temperature");
      gr2->Draw("AL");
      temperature.cd(3);
      TGraph *grd = new TGraph(n,xDif,yDif);
      grd->SetTitle("Dif temperature");
      grd->Draw("AL");
    } 
    
  
    //if ( asu2 )  {  temperature.cd(2); asu2->Draw(); } 
    //if ( dif )  {  temperature.cd(3); dif->Draw(); } 
    temperature.Write();
  }


  statistic->Write();
  deltaT.Write();

  flagTree->SetFillColor(45);
  TCanvas flagCanvas("eventFlags", "event flag value");
  flagCanvas.SetFillColor(42);
  flagCanvas.SetGrid();
  flagCanvas.SetLogx();
  flagCanvas.SetFrameFillColor(33);
  flagTree->Draw("Flag","", "hbar2");
  flagCanvas.Write();

  output->Close();
  f.Close(); 

  cout <<endl << endl << "Done" << endl;
  
return 0;
}


  

---