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

#include <map>
#include <iostream>
#include "TFile.h"
#include <TROOT.h>
#include <TRint.h>
#include <TKey.h>
#include "TTree.h"
#include "root/MTRun.hh"
#include "root/MTEvent.hh"
#include "root/MTChannel.hh"

Include dependency graph for border_event_tagging.cpp:

Go to the source code of this file.

Functions

int main (void)


Function Documentation

int main ( void   ) 

Definition at line 22 of file border_event_tagging.cpp.

00022               {  
00023   
00024   //-----------------
00025   //User's parameters
00026   
00027   //Do you want some info ? (0 = not at all, 1 = yes a little, 2 = yes a lot, 3 = yes all : for debug purpose)
00028   Int_t blabla = 1;
00029   
00030   //threshold (number of hits in a given asu at a given time (maximum=160))
00031   Int_t threshold = 80;
00032 
00033   //input data
00034   //string file = "05082011_1055_MGT.root";  
00035   string file = "06082011_0829_MGT.root";  
00036 
00037   //Micromegas chamber
00038   UInt_t chb = 1;
00039   //-----------------
00040   
00041 
00042   //our data file
00043   TFile *f = new TFile(file.c_str());
00044   if(blabla==1||blabla==2||blabla==3){
00045     cout<<endl;
00046     cout<<"reading File : "<<file<<endl;
00047   }
00048   
00049   //key iterator
00050   TIter nextkey(f->GetListOfKeys());
00051   TKey *key;
00052   UInt_t counter = 0;
00053    
00054   if(blabla==1||blabla==2||blabla==3) cout<<endl;
00055   
00056 
00057   //loop on all the keys
00058   while (key = (TKey*)nextkey()) {
00059     
00060     counter++;
00061     if(blabla==1||blabla==2||blabla==3) cout<<"key number = "<<counter<<endl;
00062     
00063     TTree *tree = (TTree*)key->ReadObj();                
00064     MTRun* run = (MTRun*)tree->GetUserInfo()->FindObject("MTRun");
00065 
00066     Int_t nEvent = tree->GetEntries();
00067     if(blabla==1||blabla==2||blabla==3) cout <<"Nevent = "<<nEvent<<endl<<endl;
00068     
00069     MTEvent *evt = new MTEvent();
00070     
00071     TBranch *branch= tree->GetBranch("MTEvent");
00072     branch->SetAddress(&evt);
00073     
00074     //at each new key, declare the event counters and indicators
00075     bool IsBorder;
00076     Int_t Asu;
00077     ULong_t AbsTime;
00078 
00079     
00080     //loop on all the events of the key
00081     for (int i=0;i<nEvent-1;i++){
00082       
00083       tree->GetEntry(i);
00084       
00085       //an event can be border only if it is not square
00086       UInt_t A;
00087       UInt_t T;
00088       if(evt->IsSquare(450,A,T,chb)) {
00089         if(blabla==2||blabla==3) cout<<"evt "<<i+1<<" square"<<endl;
00090         continue;
00091       }
00092       
00093       UInt_t nchannel = evt->GetNchannel();
00094       if(blabla==2||blabla==3) cout<<"Event "<<i+1<<" / "<<nEvent<<"; nb of channels=  "<<nchannel<<endl;
00095       
00096       //at each event declare the hits counter
00097       map<unsigned int,map<unsigned long,unsigned int> > channels_map;// key1=asu (board_id), key2=hit_time      
00098       
00099       //at each event reset the square indicator
00100       IsBorder=0;
00101       Asu = -1;
00102       AbsTime = -1;
00103 
00104       //loop on all the channels of the event
00105       for (unsigned int j=0;j<nchannel;j++){
00106 
00107         //this channel
00108         MTChannel* channel = (MTChannel*)evt->GetChannels()->UncheckedAt(j);
00109 
00110         if (channel->GetChamberId()==chb)
00111           { 
00112             //at each channel get board_id and hit_time
00113             ULong_t hit_time =  channel->GetBcIdAbs() - (channel->GetBcIdDif() - channel->GetBcIdHit() );
00114             UInt_t board_id = channel->GetBoardId();
00115             if(blabla==3) cout<<"Channel : "<<j<<"; hit on asu "<<board_id<<"; at time "<<hit_time<<endl;
00116             
00117             //get the position of this hit/channel
00118             //Int_t xpos = channel->GetX()/1000;
00119             //Int_t ypos = channel->GetY()/1000;
00120             if(blabla==3) cout <<"X= "<<channel->GetX()<<"; Y= "<<channel->GetY()<<endl;
00121             
00122             //if the hit is not on the borders of an asu go to next hit
00123             if (channel->GetX() != 0000 && channel->GetX() != 310000 && channel->GetX() != 320000 && channel->GetX() != 640000 && channel->GetX() != 650000 && channel->GetX() != 960000 && channel->GetY() != 000000 && channel->GetY() != 470000 && channel->GetY() != 480000 && channel->GetY() != 960000) continue;
00124             
00125             //upgrade the hits counter of the good board at the good time
00126             if (channels_map.find(board_id) == channels_map.end()) { // if this asu doesn't exist yet
00127               
00128               map<unsigned long, unsigned int> nbHitMap;
00129               nbHitMap[hit_time]=1;
00130               channels_map.insert(make_pair(board_id,nbHitMap));
00131             }
00132             else{ // if this asu already exist 
00133               
00134               map<unsigned long, unsigned int> &nbHitMap =  channels_map.find(board_id)->second;
00135               if (nbHitMap.find(hit_time) == nbHitMap.end()){ //if this hit_time doesn't exist yet
00136                 
00137                 nbHitMap.insert(make_pair(hit_time,1)); 
00138               }
00139               else{ //if this abs time already exist
00140                 
00141                 nbHitMap[hit_time] = nbHitMap[hit_time]++;
00142               }
00143             }  
00144           }      
00145         
00146         //Loop on each board, at each time, to identify if this event is border
00147         map<unsigned int, map<unsigned long,unsigned int> >::iterator iter1;
00148         for(iter1=channels_map.begin();iter1!=channels_map.end();iter1++){
00149           
00150           unsigned int azu=iter1->first;
00151           map<unsigned long,unsigned int> hits_per_time=iter1->second;
00152           
00153           map<unsigned long,unsigned int>::iterator iter2;
00154           for(iter2=hits_per_time.begin();iter2!=hits_per_time.end();iter2++){
00155             
00156             unsigned long time = iter2->first;
00157             unsigned int nb_hits = iter2->second;
00158             
00159             //if this is a border event
00160             if (nb_hits>=threshold) {
00161               
00162               //set square event indicator to 1
00163               IsBorder = 1;
00164               Asu = azu;
00165               AbsTime = time;
00166               
00167               if(blabla==1||blabla==2||blabla==3) cout <<" -> event: "<<i+1<<" = border in asu no "<<azu<<" at abstime "<<time<<" with n_hits = "<<nb_hits<<endl;
00168             }
00169           }
00170         }
00171       }
00172     }
00173     if(blabla==2||blabla==3) cout<<endl;
00174   }
00175   if(blabla==1||blabla==2||blabla==3) cout<<endl;
00176   return 0; 
00177 }


Generated on Mon Jan 7 13:17:09 2013 for MicromegasFramework by  doxygen 1.4.7