Ignore:
Timestamp:
01/04/15 05:11:16 (10 years ago)
Author:
Daniela Dorner
Message:
added the possibility to call macro also for 1 ganymed file
File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/Mars/fact/processing/numevents.C

    r17977 r18075  
    88#include <TSQLResult.h>
    99#include <TSQLRow.h>
     10#include <TPRegexp.h>
    1011
    1112#include "MSQLMagic.h"
    … …  
    1718using namespace std;
    1819
    19 //missing:
    20 // remove sourcekey in paths
    21 // create tables in DB
    22 
    23 int numevents(TString night="20130418", TString inpath="/daq/analysis/1/", TString table= "AnalysisResultsRunLP", Bool_t dummy=kTRUE, Bool_t pernight=kFALSE, Int_t source=0)
     20int process(MSQLMagic &serv, TString gfile, TString sfile, TString night2, TString runid, TString table, Bool_t pernight, Int_t source)
     21{
     22    TString where;
     23    TString vars;
     24    //the variables to save into the db
     25    Int_t fNumEvtsAfterCleaning = 0;
     26    Int_t fNumEvtsAfterQualCuts = 0;
     27    Int_t fNumEvtsAfterBgCuts   = 0;
     28    Float_t fNumBgEvts          = 0;
     29    Float_t fNumSigEvts         = 0;
     30    Float_t fNumExcEvts         = 0;
     31    Float_t fNumIslandsMean     = 0;
     32    Float_t fOnTimeAfterCuts    = 0;
     33
     34    if (!pernight)
     35    {
     36        //check star file status
     37        if (gSystem->AccessPathName(sfile))
     38        {
     39            cout << "ERROR - file " << sfile << " does not exist." << endl;
     40            return 2;
     41        }
     42
     43        //try to open the star file
     44        TFile star_file(sfile, "READ");
     45        if (!star_file.IsOpen())
     46        {
     47            cout << "ERROR - Could not open file " << sfile << endl;
     48            return 2;
     49        }
     50
     51        //Get the events tree
     52        TTree *star_tree = (TTree *)star_file.Get("Events");
     53        if (!star_tree)
     54        {
     55            cout << "ERROR - Could not read tree " << sfile << endl;
     56            return 2;
     57        }
     58
     59        //the number of events after cleaning
     60        fNumEvtsAfterCleaning = star_tree->GetEntries();
     61
     62        //get the mean number of islands in the _I.root file on runbasis
     63        star_tree->Draw("MImagePar.fNumIslands>>IslandHisto","","");
     64        TH1F *HistJo = (TH1F*)gDirectory->Get("IslandHisto");
     65        if (!HistJo)
     66        {
     67            cout << "ERROR - Reading of IslandHisto failed " << HistJo << endl;
     68            return 2;
     69        }
     70        fNumIslandsMean = HistJo->GetMean();
     71    }
     72
     73    //check ganymed file status
     74    if (gSystem->AccessPathName(gfile))
     75    {
     76        cout << "ERROR - file " << gfile << " does not exist." << endl;
     77        return 2;
     78    }
     79
     80    //try to open the ganymed file
     81    TFile ganymed_file(gfile, "READ");
     82    if (!ganymed_file.IsOpen())
     83    {
     84        cout << "ERROR - Could not open file " << gfile << endl;
     85        return 2;
     86    }
     87
     88    //Get the status display contents
     89    MStatusArray arr;
     90    if (arr.Read()<=0)
     91    {
     92        cout << "ERROR - could not read MStatusArray, file " << gfile << endl;
     93        return 2;
     94    }
     95
     96    TH1D *vstime = (TH1D*)arr.FindObjectInCanvas("OnTime",  "TH1D", "OnTime");
     97    if (!vstime)
     98    {
     99        cout << "WARNING - Reading of Theta failed." << endl;
     100        return 2;
     101    }
     102    fOnTimeAfterCuts = vstime->Integral();
     103
     104    //Get number of events after quality cuts
     105    //from the number of entries in a histogram
     106    TH1F *precut =
     107        (TH1F*)arr.FindCanvas("PreCut")->FindObject("SizeSame");
     108    if (!precut)
     109    {
     110        cout << "WARNING - Reading of PreCut canvas failed." << endl;
     111        return 2;
     112    }
     113    fNumEvtsAfterQualCuts = precut->GetEntries();
     114
     115    //Get number of events after background rejection cuts
     116    //from the number of entries in a histogram
     117    TH1F *postcut =
     118        (TH1F*)arr.FindCanvas("PostCut")->FindObject("SizeSame");
     119    if (!postcut)
     120    {
     121        cout << "WARNING - Reading of PostCut canvas failed." << endl;
     122        return 2;
     123    }
     124    fNumEvtsAfterBgCuts = postcut->GetEntries();
     125
     126    //Get signal, bg and excess events
     127    MHThetaSq *halphaon = (MHThetaSq*)arr.FindObjectInCanvas("MHThetaSq", "Hist");
     128    if (!halphaon)
     129    {
     130        cout << "WARNING - Reading of MHThetaSq failed, file " << gfile << endl;
     131        return 2;
     132    }
     133    const MAlphaFitter &fit = halphaon->GetAlphaFitter();
     134    if (!&fit)
     135    {
     136        cout << "WARNING - Reading of MAlphaFitter failed, file " << gfile << endl;
     137        return 2;
     138    }
     139
     140    //check if the excess events number is small and close to 0
     141    //due to the computer precision error when subtracting sig and bg events
     142    fNumBgEvts  = fit.GetEventsBackground();
     143    fNumSigEvts = fit.GetEventsSignal();
     144    fNumExcEvts = fit.GetEventsExcess();
     145    fNumExcEvts = fabs(fNumExcEvts)<1e-5 ? 0 : fNumExcEvts;
     146
     147    /*
     148     cout<< fNumEvtsAfterCleaning<< "    "
     149     << fNumEvtsAfterQualCuts<< "    "
     150     << fNumEvtsAfterBgCuts  << "    "
     151     << fNumBgEvts           << "    "
     152     << fNumSigEvts          << "    "
     153     << fNumExcEvts          << endl;
     154     */
     155
     156    //inserting or updating the information in the database
     157    if (!pernight)
     158        vars = Form("fRunID=%s, fNight=%s,", runid.Data(), night2.Data());
     159    else
     160        vars = Form("fNight=%s, fSourceKey=%d, ", night2.Data(), source);
     161    vars += Form(" fNumEvtsAfterQualCuts=%d, "
     162                 " fNumEvtsAfterBgCuts=%d, "
     163                 " fNumBgEvts=%6.1f, "
     164                 " fNumSigEvts=%6.1f, "
     165                 " fNumExcEvts=%6.1f, "
     166                 " fOnTimeAfterCuts=%7.2f ",
     167                 fNumEvtsAfterQualCuts,
     168                 fNumEvtsAfterBgCuts,
     169                 fNumBgEvts,
     170                 fNumSigEvts,
     171                 fNumExcEvts,
     172                 fOnTimeAfterCuts
     173                );
     174    if (!pernight)
     175        vars += Form(", fNumIslandsMean=%6.2f, fNumEvtsAfterCleaning=%d ",
     176                     fNumIslandsMean,
     177                     fNumEvtsAfterCleaning
     178                    );
     179
     180
     181    if (pernight)
     182        where = Form("fNight=%s AND fSourceKey=%d", night2.Data(), source);
     183    else
     184        where = Form("fRunID=%s AND fNight=%s", runid.Data(), night2.Data());
     185
     186    cout << "vars: " << vars << endl;
     187    cout << "where: " << where << endl;
     188
     189    if (!serv.InsertUpdate(table, vars, where))
     190    {
     191        cout << "ERROR - insert/update to DB failed  (vars:" << vars << " where: " << where << ")." << endl;
     192        return 2;
     193    }
     194    return 1;
     195}
     196
     197
     198int numevents(Int_t night, TString inpath="/loc_data/analysis/", TString table= "AnalysisResultsRunLP", Bool_t dummy=kTRUE, Bool_t pernight=kFALSE, Int_t source=0)
    24199{
    25200
    … …  
    41216    //and only data runs
    42217    query+=" WHERE fRunTypeKEY=1 ";
    43     query+=" AND fNight="+night;
     218    query+=Form(" AND fNight=%d",night);
    44219    //query+=" AND fRunID=128 ";
    45220    if (pernight)
    … …  
    65240    TString ganymed_fname;
    66241    TString star_fname;
    67     TString where;
    68     TString vars;
    69 
    70     //the variables to save into the db
    71     Int_t fNumEvtsAfterCleaning = 0;
    72     Int_t fNumEvtsAfterQualCuts = 0;
    73     Int_t fNumEvtsAfterBgCuts   = 0;
    74     Float_t fNumBgEvts          = 0;
    75     Float_t fNumSigEvts         = 0;
    76     Float_t fNumExcEvts         = 0;
    77     Float_t fNumIslandsMean     = 0;
    78     Float_t fOnTimeAfterCuts    = 0;
     242    Int_t rc=0;
    79243
    80244    //loop over the data files
    81245    TSQLRow *row=0;
    82246
    83     /*
    84247    //cout the data names
    85248    cout<<"#evts_after_cleaning" << "    "
    … …  
    89252        <<"evts_from_on_region"  << "    "
    90253        <<"calculated_excess_evts_in_file" << "    " << endl;
    91     */
    92254
    93255    while ((row=res->Next()))
    … …  
    115277            star_fname=inpath+"/star/"+year+"/"+month+"/"+day+"/"+run+"_I.root";
    116278            cout << "sf => " << star_fname << endl;
    117 
    118             //check star file status
    119             if (gSystem->AccessPathName(star_fname))
    120             {
    121                 cout << "ERROR - file " << star_fname << " does not exist." << endl;
    122                 continue;
    123             }
    124 
    125             //try to open the star file
    126             TFile star_file(star_fname, "READ");
    127             if (!star_file.IsOpen())
    128             {
    129                 cout << "ERROR - Could not open file " << star_fname << endl;
    130                 continue;
    131             }
    132 
    133             //Get the events tree
    134             TTree *star_tree = (TTree *)star_file.Get("Events");
    135             if (!star_tree)
    136             {
    137                 cout << "ERROR - Could not read tree " << star_fname << endl;
    138                 continue;
    139             }
    140 
    141             //the number of events after cleaning
    142             fNumEvtsAfterCleaning = star_tree->GetEntries();
    143 
    144             //get the mean number of islands in the _I.root file on runbasis
    145             star_tree->Draw("MImagePar.fNumIslands>>IslandHisto","","");
    146             TH1F *HistJo = (TH1F*)gDirectory->Get("IslandHisto");
    147             if (!HistJo)
    148             {
    149                 cout << "ERROR - Reading of IslandHisto failed " << HistJo << endl;
    150                 continue;
    151             }
    152             fNumIslandsMean = HistJo->GetMean();
    153         }
    154 
    155         //check ganymed file status
    156         if (gSystem->AccessPathName(ganymed_fname))
    157         {
    158             cout << "ERROR - file " << ganymed_fname << " does not exist." << endl;
    159             continue;
    160         }
    161 
    162         //try to open the ganymed file
    163         TFile ganymed_file(ganymed_fname, "READ");
    164         if (!ganymed_file.IsOpen())
    165         {
    166             cout << "ERROR - Could not open file " << ganymed_fname << endl;
    167             continue;
    168         }
    169 
    170         //Get the status display contents
    171         MStatusArray arr;
    172         if (arr.Read()<=0)
    173         {
    174             cout << "ERROR - could not read MStatusArray, file " << ganymed_fname << endl;
    175             continue;
    176         }
    177 
    178         TH1D *vstime = (TH1D*)arr.FindObjectInCanvas("OnTime",  "TH1D", "OnTime");
    179         if (!vstime)
    180         {
    181             cout << "WARNING - Reading of Theta failed." << endl;
    182             continue;
    183         }
    184         fOnTimeAfterCuts = vstime->Integral();
    185 
    186         //Get number of events after quality cuts
    187         //from the number of entries in a histogram
    188         TH1F *precut =
    189                 (TH1F*)arr.FindCanvas("PreCut")->FindObject("SizeSame");
    190         if (!precut)
    191         {
    192             cout << "WARNING - Reading of PreCut canvas failed." << endl;
    193             continue;
    194         }
    195         fNumEvtsAfterQualCuts = precut->GetEntries();
    196 
    197         //Get number of events after background rejection cuts
    198         //from the number of entries in a histogram
    199         TH1F *postcut =
    200                 (TH1F*)arr.FindCanvas("PostCut")->FindObject("SizeSame");
    201         if (!postcut)
    202         {
    203             cout << "WARNING - Reading of PostCut canvas failed." << endl;
    204             continue;
    205         }
    206         fNumEvtsAfterBgCuts = postcut->GetEntries();
    207 
    208         //Get signal, bg and excess events
    209         MHThetaSq *halphaon = (MHThetaSq*)arr.FindObjectInCanvas("MHThetaSq", "Hist");
    210         if (!halphaon)
    211         {
    212             cout << "WARNING - Reading of MHThetaSq failed, file " << ganymed_fname << endl;
    213             continue;
    214         }
    215         const MAlphaFitter &fit = halphaon->GetAlphaFitter();
    216         if (!&fit)
    217         {
    218             cout << "WARNING - Reading of MAlphaFitter failed, file " << ganymed_fname << endl;
    219             continue;
    220         }
    221 
    222         //check if the excess events number is small and close to 0
    223         //due to the computer precision error when subtracting sig and bg events
    224         fNumBgEvts  = fit.GetEventsBackground();
    225         fNumSigEvts = fit.GetEventsSignal();
    226         fNumExcEvts = fit.GetEventsExcess();
    227         fNumExcEvts = fabs(fNumExcEvts)<1e-5 ? 0 : fNumExcEvts;
    228 
    229         /*
    230         cout<< fNumEvtsAfterCleaning<< "    "
    231             << fNumEvtsAfterQualCuts<< "    "
    232             << fNumEvtsAfterBgCuts  << "    "
    233             << fNumBgEvts           << "    "
    234             << fNumSigEvts          << "    "
    235             << fNumExcEvts          << endl;
    236         */
    237 
    238         //inserting or updating the information in the database
    239         if (!pernight)
    240             vars = Form("fRunID=%s, fNight=%s,", runid.Data(), night2.Data());
    241         else
    242             vars = Form("fNight=%s, fSourceKey=%d, ", night2.Data(), source);
    243         vars += Form(" fNumEvtsAfterQualCuts=%d, "
    244                     " fNumEvtsAfterBgCuts=%d, "
    245                     " fNumBgEvts=%6.1f, "
    246                     " fNumSigEvts=%6.1f, "
    247                     " fNumExcEvts=%6.1f, "
    248                     " fOnTimeAfterCuts=%7.2f ",
    249                     fNumEvtsAfterQualCuts,
    250                     fNumEvtsAfterBgCuts,
    251                     fNumBgEvts,
    252                     fNumSigEvts,
    253                     fNumExcEvts,
    254                     fOnTimeAfterCuts
    255                    );
    256         if (!pernight)
    257             vars += Form(", fNumIslandsMean=%6.2f, fNumEvtsAfterCleaning=%d ",
    258                          fNumIslandsMean,
    259                          fNumEvtsAfterCleaning
    260                         );
    261                    
    262 
    263         if (pernight)
    264             where = Form("fNight=%s AND fSourceKey=%d", night2.Data(), source);
    265         else
    266             where = Form("fRunID=%s AND fNight=%s", runid.Data(), night2.Data());
    267 
    268         cout << "vars: " << vars << endl;
    269         cout << "where: " << where << endl;
    270 
    271         if (!serv.InsertUpdate(table, vars, where))
    272         {
    273             cout << "ERROR - insert/update to DB failed  (vars:" << vars << " where: " << where << ")." << endl;
    274             return 2;
    275         }
    276 
    277     }
     279        }
     280        rc=process(serv, ganymed_fname, star_fname, night2, runid, table, pernight, source);
     281        if (rc>1)
     282            return rc;
     283
     284    }
     285
    278286
    279287    return 1;
    280288}
     289
     290int numevents(TString ganymedfile, TString starfile, TString table, Bool_t dummy=kTRUE)
     291{
     292
     293    //connect to mysql server
     294    MSQLMagic serv("sql.rc");
     295    if (!serv.IsConnected())
     296    {
     297        cout << "ERROR - Connection to database failed." << endl;
     298        return 0;
     299    }
     300    serv.SetIsDummy(dummy);
     301    //get night and run from file name
     302    TString basename=gSystem->BaseName(ganymedfile);
     303    TPRegexp regexp("20[0-9][0-9][0-2][0-9][0-3][0-9]_[0-9][0-9][0-9]");
     304    TString run = basename(regexp);
     305    TPRegexp regexp2("20[0-9][0-9][0-2][0-9][0-3][0-9]");
     306    TPRegexp regexp3("_[0-9][0-9][0-9]");
     307    TPRegexp regexp4("[0-9][0-9][0-9]");
     308    TString night = run(regexp2);
     309    TString runnum = run(regexp3);
     310    TString runid = runnum(regexp4);
     311
     312    return process(serv, ganymedfile, starfile, night, runid, table, false, 0);
     313}
     314
Note: See TracChangeset for help on using the changeset viewer.