/data3/calcul/jacquem/working_dir/Micromegas/micromegasFrameWork/src/analyse/root/buildCaloHits.cpp File Reference

#include <TRandom.h>
#include <TTree.h>
#include <TFile.h>
#include <TStyle.h>
#include <TKey.h>
#include <Riostream.h>
#include <sstream>
#include <string>
#include "mysql/Mysql.hh"
#include "MicroException.hh"
#include <iostream>
#include <fstream>
#include "root/CaloHit.hh"
#include "TH2I.h"
#include "root/CaloRun.hh"
#include "root/CaloEvent.hh"
#include "root/MTRun.hh"
#include "root/MTEvent.hh"
#include "root/MTChannel.hh"
#include "root/MTTrack.hh"

Include dependency graph for buildCaloHits.cpp:

Go to the source code of this file.

Functions

void selectAndCopyEvent (const string &inputFileName, const string &outputFileName, int deltaTMin, int deltaTMax)
void buildCaloEvent (MTEvent &evt, CaloEvent &, int deltaTMin, int deltaTMax, float offset_x, float offset_y)
void usage (void)
int main (int argc, char **argv)


Function Documentation

void selectAndCopyEvent ( const string &  inputFileName,
const string &  outputFileName,
int  deltaTMin,
int  deltaTMax 
)

Definition at line 156 of file buildCaloHits.cpp.

00157 {
00158 
00159  
00160   TFile readFile(inputFileName.c_str(),"READONLY");
00161   TFile caloFile(outputFileName.c_str(), "UPDATE","Micromegas calorimeter data");
00162 
00163   ui16 runId = 0;
00164 
00165 
00166   TIter nextkey(readFile.GetListOfKeys());
00167   TKey *key = NULL;
00168   CaloRun run;
00169 
00170   key= (TKey*)nextkey();
00171   
00172 
00173   //Get tree to read
00174   TTree *readTree = (TTree*)key->ReadObj();
00175   MTEvent *readEvt =  new MTEvent();
00176   TBranch *branch= readTree->GetBranch("MTEvent");
00177   branch->SetAddress(&readEvt);
00178 
00179 
00180   //Get RunId
00181   MTRun * readRun  = (MTRun*)readTree->GetUserInfo()->FindObject("MTRun");
00182   int runid = readRun->GetRunId();
00183   if (runid==0)
00184     {
00185       //cout<<"RunId not known, please enter: ";
00186       //cin>>runid;
00187       //cout<<endl;
00188       //cout<<endl;
00189       runid = 10227;
00190     }
00191 
00192   //Get telescope/prototype offset from database
00193   Mysql mysql;
00194   float offset_x = mysql.getRunOffsetX(runid);
00195   float offset_y = mysql.getRunOffsetY(runid);
00196 
00197   cout<<"OFFSETS = "<<offset_x<<" "<<offset_y<<endl;
00198 
00199 
00200 
00201   // Create a new TTree which will contain  filtered events
00202   stringstream label;
00203   label<< inputFileName << "_" << runId ; 
00204   TTree tree(label.str().c_str(),"Calo tree");
00205   cout << "work on :" << label.str() << endl;
00206   tree.SetMaxTreeSize(4000000000); // 4 Mb
00207 
00208   CaloEvent *evt = new CaloEvent();
00209   // create a branch with run 
00210   tree.Branch("CaloEvent",&evt,16000,2);
00211   run.SetRunId(runId++);
00212 
00213 
00214   // Loop over event in current input TTree
00215   int nbEntries = readTree->GetEntries();
00216   ui16 evtNum = 0;
00217   cout<<"Tree entries = "<<nbEntries<<endl;
00218 
00219   MTTrack * trk = NULL;
00220 
00221 
00222   //for ( evtNum = 0; evtNum < 100 ; evtNum++)
00223   for ( evtNum = 0; evtNum < nbEntries ; evtNum++)
00224     {
00225       cout << "Done for entry " << evtNum + 1 << " / " << nbEntries << "\r" << flush;
00226 
00227       // Get input Event as readEvt
00228       readTree->GetEntry(evtNum);
00229 
00230       bool valid_time = readEvt->GetTimeInfo();
00231       bool square_m2 = readEvt->GetSquareEvt();
00232 
00233       if (valid_time && !square_m2)
00234         {
00235 
00236           trk = readEvt->GetTrack();
00237 
00238           if (trk != NULL)
00239             {
00240               buildCaloEvent(*readEvt,*evt,deltaTMin,deltaTMax,offset_x,offset_y);
00241               tree.Fill();
00242               evt->Clear();
00243             }
00244         }
00245     }
00246 
00247 
00248   tree.GetUserInfo()->Add(&run); 
00249   TFile *curFile =tree.GetCurrentFile() ;
00250   //curFile->Write("", TObject::kOverwrite);
00251   tree.Write();
00252 
00253   delete readEvt;
00254   delete evt;
00255   
00256 
00257   caloFile.Close();
00258   readFile.Close();
00259   // delete run;
00260 }

void buildCaloEvent ( MTEvent evt,
CaloEvent ,
int  deltaTMin,
int  deltaTMax,
float  offset_x,
float  offset_y 
)

Definition at line 219 of file buildCaliceCaloHits.cpp.

00220 {
00221 
00222     /*m2 prototype time cut*/
00223     int dt_min = 5;
00224     int dt_max = 10;
00225 
00226     if (deltaTMin != -1) { deltaTMin = dt_min; }
00227     if (deltaTMax != -1) { deltaTMax = dt_max; }
00228 
00229 
00230     //Distance from first telescope chamber to m2 prototype
00231     float dist_z = 156.;
00232 
00233 
00234     //Get track
00235     MTTrack * trk = evt.GetTrack();
00236     
00237     int trk_quality = 0;
00238     if ( trk != NULL )
00239     {
00240       if (trk->GetFromPad())
00241       {
00242         trk_quality = 10;
00243       }
00244 
00245       if (trk->GetFromStrip())
00246       {
00247         trk_quality = 20;
00248       }
00249 
00250       if (trk->GetStraight())
00251       {
00252         trk_quality++;
00253 
00254         if (trk->GetSingleHit())
00255         {
00256          trk_quality++;
00257         }
00258       }
00259 
00260 
00261       caloEvent.SetTrkQuality(trk_quality);
00262 
00263       float extrapol_x = offset_x;
00264       float extrapol_y = offset_y;
00265 
00266       //Extrapolate at prototype
00267 
00268       //straight tracks
00269       if (trk_quality%10 != 0)
00270       {
00271       extrapol_x += trk->GetBx();
00272       extrapol_y += trk->GetBy();       
00273       }
00274 
00275       //inclined tracks
00276       else
00277       {
00278       extrapol_x += trk->GetAx() * dist_z + trk->GetBx();
00279       extrapol_y += trk->GetAy() * dist_z + trk->GetBy();
00280       }
00281 
00282 
00283       caloEvent.SetXExtrapol(extrapol_x);
00284       caloEvent.SetYExtrapol(extrapol_y);
00285   }//  end if track exist
00286 
00287 
00288   
00289 
00290     //loop over hits
00291 
00292     int nchannel = evt.GetNchannel();
00293     int nhit_m2_cut = 0;
00294 
00295     for (int j=0;j<nchannel;j++)
00296       {
00297         MTChannel* channel = (MTChannel*)evt.GetChannels()->UncheckedAt(j);
00298           
00299         //look in m2 prototype
00300         if (channel->GetChamberId() == 10)
00301           {
00302             Long64_t dt = channel->GetBcId_Dif() - channel->GetBcId_Hit();
00303               
00304             int zpos = (int)channel->GetZ() / 10000;
00305 
00306             if (dt>=dt_min && dt<=dt_max)
00307               {
00308                 int digit = channel->GetDigitalValue();
00309                 int xpos = (int)channel->GetX() / 10000;
00310                 int ypos = (int)channel->GetY() / 10000;
00311                   
00312                 int t = (int)dt;
00313                   
00314                 CaloHit hit;
00315 
00316                 hit.SetX(xpos);
00317                 hit.SetY(ypos);           
00318                 hit.SetZ(zpos);
00319                 hit.SetDeltaT(t);
00320                 hit.SetThreshold(digit);
00321 
00322                 caloEvent.AddHit(hit);  
00323    }
00324   }
00325  }
00326 }

void usage ( void   ) 

Definition at line 32 of file buildCaloHits.cpp.

00033 {
00034     cout  << "usage: buildCaloHits -t timestampMin timestampMax outputfileName [deltaTMin deltaTMax]" << endl;
00035     cout  << "     : buildCaloHits -f filename outputfileName [deltaTMin deltaTMax]" << endl;
00036     cout << "1312577280 1312590660" << endl;
00037 }

int main ( int  argc,
char **  argv 
)

Definition at line 39 of file buildCaloHits.cpp.

00039                                 {
00040 
00041   if ( argc < 4 || argc > 7 )
00042   {
00043     usage();
00044     return -1;
00045   }
00046 
00047 
00048   int deltaTMin = -1;
00049   int deltaTMax = -1;
00050 
00051 
00052 
00053 
00054 
00055   // Get inputFile name in database by timestamp limit
00056   if ( strcmp (argv[1] ,"-t" )== 0 )
00057   {
00058     if ( argc < 5 ) { usage(); return -1; } 
00059 
00060     if ( argc == 7 ) 
00061     {
00062       deltaTMin = atoi ( argv[5] );
00063       deltaTMax = atoi ( argv[6] );
00064     }
00065     
00066     string outputFileName = argv[4];
00067     // RECREATE Root file to clear its contains if it alread exist
00068     TFile caloFile(outputFileName.c_str(), "RECREATE","Micromegas calorimeter data");
00069     caloFile.Write();
00070     caloFile.Close();
00071 
00072     Mysql sql;
00073     string start = argv[2];
00074     string end = argv[3];
00075     
00076     try 
00077     {
00078       // Get list of TFile to filter from Database
00079       const vector<string>& ids = sql.getMergedFilesByTimestamp(start,end);
00080       vector<string>::const_iterator cii;
00081 
00082       // Loop over input TFile
00083       for(cii=ids.begin(); cii!=ids.end(); cii++)
00084       {
00085 
00086         string inputFileName = *cii;
00087         selectAndCopyEvent(inputFileName,outputFileName,deltaTMin,deltaTMax);
00088       }
00089       delete &ids;
00090     }
00091     catch (MicroException e) { }
00092   }
00093 
00094 
00095 
00096 
00097 
00098 
00099 
00100   // Get input file name in a input file
00101   else if ( strcmp (argv[1] ,"-f") == 0 )
00102   {
00103      if ( argc < 4 ) { usage(); return -1; } 
00104      string outputFileName = argv[3];
00105 
00106     if ( argc == 6 ) 
00107     {
00108       deltaTMin = atoi ( argv[4] );
00109       deltaTMax = atoi ( argv[5] );
00110     }
00111 
00112     
00113 
00114 
00115 
00116     // RECREATE Root file to clear its contains if it alread exist
00117     TFile caloFile(outputFileName.c_str(), "RECREATE","Micromegas calorimeter data");
00118     caloFile.Write();
00119     caloFile.Close();
00120 
00121     ifstream myReadFile;
00122     myReadFile.open(argv[2]);
00123     string inputRootFile;
00124     if (myReadFile.is_open()) {
00125       while (!myReadFile.eof()) {
00126         getline(myReadFile,inputRootFile);
00127         if ( inputRootFile.size() > 0 ) 
00128         {
00129           cout<<"Read ROOT file : "<<inputRootFile<<endl;
00130           selectAndCopyEvent(inputRootFile,outputFileName,deltaTMin,deltaTMax);
00131         }
00132       }
00133     }
00134     myReadFile.close();
00135     return 0;
00136   }
00137   else 
00138   {
00139     usage();
00140     return -1;
00141   }
00142   return 0;
00143 }


Generated on Mon Jun 11 16:57:15 2012 for MicromegasFramework by  doxygen 1.4.7