/data3/calcul/jacquem/working_dir/Micromegas/micromegasFrameWork/src/analyse/root/Tag_ramfull.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) $ */
#include "tools/Log.hh"
#include "tools/Toolbox.hh"

#include "root/MTRun.hh"
#include "root/MTChannel.hh"
#include "root/MTEvent.hh"
#include "root/MTTrack.hh"
#include "MicroException.hh"

#include "mysql/Mysql.hh"


#include "mysql/Mysql.hh"

#include "TFile.h"
#include "TTree.h"
#include "TKey.h"
#include "TBranch.h"
#include "TRandom.h"
#include "TProfile.h"
#include "TH2I.h"
#include "TF1.h"


#include <sstream>
#include <sys/types.h>
#include <unistd.h>
#include <stdio.h>
#include <stdlib.h>


int errsv = errno;


using namespace std;

/*
 Modifying a existing TTree is often not possible because the tree is in a read-only
 file and you do not have permission to save the modified tree .
 Even if you do have the permission, you risk losing the
 original tree with an unsuccessful attempt to save  the modification.
*/

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

   if ( argc != 3  ) {
       FILE_LOG(logERROR)  << "usage:" << endl;
       FILE_LOG(logERROR)  << ":Tag_ramfull fileToProcess threshold" << endl;
       exit(1);
   }
  string inputfile = argv[1] ;

  UInt_t threshold = atoi(argv[2]);

  // build temporaly fill path and name
  pid_t pid = getpid(); 
  std::ostringstream oss;
  oss << "/tmp/" << pid << ".root";
  string outputfile = oss.str();

  cout<<"Read   "<<inputfile<<endl;



  // Open Root TFile in read only mode
  TFile infile(inputfile.c_str(),"READ");

  // Open Root TFile in write mode
  TFile oufile(outputfile.c_str(),"RECREATE");
  // Change current directory to "this" directory
  // Mandatory to use updateTree->GetCurrentFile() later
  oufile.cd();



  // Loop over original TTrees
  TIter nextkey(infile.GetListOfKeys());
  MTEvent *updateEvt = new MTEvent();

  TKey* key = NULL;//(TKey*)nextkey();
  TTree *updateTree = NULL;
  while (key = (TKey*)nextkey())
  //key = (TKey*)nextkey();
  {
    TTree *tree = (TTree*)key->ReadObj();

    // Create new TTree which will contain processed MTRun and MTEvent
    string treeName = tree->GetName();
    treeName = treeName + "_preprocessed";
    updateTree = new TTree(treeName.c_str(),"MicroMegas event");
    updateTree->Branch("MTEvent",&updateEvt,16000,2);

    MTRun * run  = (MTRun*)tree->GetUserInfo()->FindObject("MTRun");

    // Set Processed flag to MTRun
    run->SetProcessed(true);

    // Save MTRun information to new processed TTree
    updateTree->GetUserInfo()->Add(run);
    MTEvent *evt =  new MTEvent();
    TBranch *branch= tree->GetBranch("MTEvent");
    branch->SetAddress(&evt);


    //MTChannel* channel =NULL;
    int nbEntries = tree->GetEntries();
    MTEvent* prev_evt = NULL;

    cout<<endl;
    cout<<"NUMBER OF EVENTS = "<<nbEntries<<endl;
    cout<<endl;


    for ( int evtNum = 0; evtNum < nbEntries ; evtNum++)
    {
      map<UInt_t, UInt_t> dtMap; // key is ( BcIdDif() - BcIdHit. Value is number of entry
      tree->GetEntry(evtNum);
//cout << "tree->GetEntry("<< evtNum <<")" << endl;

      // Copy original MTEvent to new processed MTEvent
      *updateEvt = *evt;

      MTChannel* channel = NULL;
      int nbChannel = evt->GetNchannel();
      // Compute dt
      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;

        // Add entry for this 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]++;
          //cout << "Event [" << evtNum << "] BcIdDif[" << t2 << "], BcIdHit[" << t3 << "], dt[" << dt << "], count " << dtMap[dt] << endl;
        }
      } // end for each channel


      std::vector<UInt_t>& dtValid = updateEvt->GetDtCandidate();
      dtValid.clear();
      for ( map<UInt_t, UInt_t>::const_iterator iter= dtMap.begin(); iter != dtMap.end(); iter++)
      {
        UInt_t dt = (*iter).first;
        UInt_t value = (*iter).second;
      
        if ( value > threshold )
        {
          dtValid.push_back(dt);
//cout << "cout <<  dtValid.push_back( "<< dt << ")" << endl;
        }
      }
      sort(dtValid.begin(),dtValid.end());
      //Write MTEvent into Branch
      cout  << "Write Event[" << evtNum << "] \r" << flush;
      updateTree->Fill();
    }
    // free heap memory
    //delete evt;
    // Write Root Object to TFile
    TFile * curFile = updateTree->GetCurrentFile();
    updateTree->Write("", TObject::kOverwrite);
    updateTree->SetDirectory(0);
  }




  infile.Close();
  oufile.Close();



  // Replace original root file by the new processed one
  string command = "cp -f ";
  command = command + outputfile + " "  + inputfile ;
  int result = system(command.c_str());
  if ( result == -1 ) { cout <<  strerror(errno) << endl; }     
  command = "rm -f " + outputfile;
  result = system(command.c_str());
  if ( result == -1 ) { cout <<  strerror(errno) << endl; }     

  cout << endl << "Done."  << endl;


return 0;
}

---