Ignore:
Timestamp:
11/21/05 11:09:12 (19 years ago)
Author:
tbretz
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars/mranforest
Files:
8 edited

Legend:

Unmodified
Added
Removed

  • trunk/MagicSoft/Mars/mranforest/MHRanForestGini.cc

    r2296 r7413  
    4343#include "MRanTree.h"
    4444#include "MRanForest.h"
     45#include "MDataArray.h"
    4546
    4647#include "MLog.h"
    … …  
    5758//
    5859MHRanForestGini::MHRanForestGini(Int_t nbins, const char *name, const char *title)
     60    : fRules(0.01, 0.01, 0.99, 0.99)
    5961{
    6062    //
    … …  
    6466    fTitle = title ? title : "Measure of importance of Random Forest-input parameters";
    6567
    66     fGraphGini = new TGraph;
    67     fGraphGini->SetTitle("Importance of RF-input parameters measured by mean Gini decrease");
    68     fGraphGini->SetMarkerStyle(kFullDotSmall);
    69 }
    70 
    71 // --------------------------------------------------------------------------
    72 //
    73 // Delete the histograms.
    74 //
    75 MHRanForestGini::~MHRanForestGini()
    76 {
    77     delete fGraphGini;
     68    fGraphGini.SetNameTitle("Gini", "Importance of RF-input parameters measured by mean Gini decrease");
     69    fGraphGini.SetMarkerStyle(kFullDotMedium);
     70
     71    fGraphError.SetNameTitle("ResErr", "Resolution/Error versus train step");
     72    fGraphError.SetMarkerStyle(kFullDotMedium);
     73
     74    fGraphNodes.SetNameTitle("Nodes", "Number of nodes versus train step");
     75    fGraphNodes.SetMarkerStyle(kFullDotMedium);
     76
     77    fRules.SetTextAlign(13);
     78    fRules.SetTextSize(0.05);
    7879}
    7980
    … …  
    104105Bool_t MHRanForestGini::Fill(const MParContainer *par, const Stat_t w)
    105106{
     107    MRanTree *t = fRanForest->GetCurTree();
     108
    106109    for (Int_t i=0;i<fRanForest->GetNumDim();i++)
    107         fGini[i]+=fRanForest->GetCurTree()->GetGiniDec(i);
     110        fGini[i] += t->GetGiniDec(i);
     111
     112    fGraphError.SetPoint(fGraphError.GetN(), GetNumExecutions(), t->GetError());
     113    fGraphNodes.SetPoint(fGraphError.GetN(), GetNumExecutions(), t->GetNumEndNodes());
    108114
    109115    return kTRUE;
    … …  
    117123    const Int_t n = fGini.GetSize();
    118124
    119     fGraphGini->Set(n);
    120 
    121     Stat_t max=0;
    122     Stat_t min=0;
     125    fGraphGini.Set(n);
     126
    123127    for (Int_t i=0; i<n; i++)
    124128    {
    125         fGini[i] /= fRanForest->GetNumTrees();
    126         fGini[i] /= fRanForest->GetNumData();
    127 
    128         const Stat_t ip = i+1;
    129         const Stat_t ig = fGini[i];
    130 
    131         if (ig>max) max=ig;
    132         if (ig<min) min=ig;
    133 
    134         fGraphGini->SetPoint(i,ip,ig);
    135     }
    136 
    137     // This is used in root>3.04/? so that SetMaximum/Minimum can succeed
    138     fGraphGini->GetHistogram();
    139 
    140     fGraphGini->SetMaximum(1.05*max);
    141     fGraphGini->SetMinimum(0.95*min);
     129        fGini[i] /= fRanForest->GetNumTrees()*fRanForest->GetNumData();
     130        fGraphGini.SetPoint(i, i+1, fGini[i]);
     131    }
     132
     133    fRules.AddText("");
     134    const MDataArray &arr = *fRanForest->GetRules();
     135    int i;
     136    for (i=0; i<arr.GetNumEntries(); i++)
     137    {
     138        TString s;
     139        s += i+1;
     140        s += ") ";
     141        s += arr.GetRule(i);
     142        fRules.AddText(s);
     143    }
     144    for (; i<20; i++)
     145        fRules.AddText("");
    142146
    143147    return kTRUE;
    … …  
    150154void MHRanForestGini::Draw(Option_t *)
    151155{
    152     if (fGraphGini->GetN()==0)
     156    if (fGraphGini.GetN()==0)
    153157        return;
    154158
    … …  
    158162    AppendPad("");
    159163
    160     fGraphGini->Draw("ALP");
    161     pad->Modified();
    162     pad->Update();
    163 
    164     TH1 *h = fGraphGini->GetHistogram();
    165     if (!h)
    166         return;
    167 
    168     //h->GetXaxis()->SetRangeUser(0, fGini.GetSize()+1);
    169     h->SetXTitle("No.of RF-input parameter");
    170     h->SetYTitle("Mean decrease in Gini-index [a.u.]");
    171 
    172     pad->Modified();
    173     pad->Update();
    174 }
     164    pad->Divide(2,2);
     165
     166    pad->cd(1);
     167    gPad->SetBorderMode(0);
     168    gPad->SetGridx();
     169    gPad->SetGridy();
     170    fGraphGini.Draw("ALP");
     171
     172    TH1 *h = fGraphGini.GetHistogram();
     173    if (h)
     174    {
     175        h->SetXTitle("No.of RF-input parameter");
     176        h->SetYTitle("Mean decrease in Gini-index [au]");
     177        h->GetXaxis()->SetNdivisions(10);
     178    }
     179
     180    pad->cd(2);
     181    gPad->SetBorderMode(0);
     182    fGraphError.Draw("ALP");
     183    h = fGraphError.GetHistogram();
     184    if (h)
     185    {
     186        h->SetXTitle("Train step/Tree number");
     187        h->SetYTitle("Error/Resolution");
     188    }
     189
     190    pad->cd(3);
     191    gPad->SetBorderMode(0);
     192    fGraphNodes.Draw("ALP");
     193    h = fGraphNodes.GetHistogram();
     194    if (h)
     195    {
     196        h->SetXTitle("Train step/Tree number");
     197        h->SetYTitle("Number of end nodes");
     198    }
     199
     200    pad->cd(4);
     201    gPad->SetBorderMode(0);
     202    fRules.Draw();
     203}

  • trunk/MagicSoft/Mars/mranforest/MHRanForestGini.h

    r2173 r7413  
    99#include <TArrayF.h>
    1010#endif
     11#ifndef ROOT_TGraph
     12#include <TGraph.h>
     13#endif
     14#ifndef ROOT_TPaveText
     15#include <TPaveText.h>
     16#endif
    1117
    12 class TH1D;
    13 class TGraph;
    1418class MParList;
    1519class MRanForest;
    … …  
    2226
    2327    TArrayF fGini;           //!
    24     TGraph *fGraphGini;      //->
     28
     29    TGraph  fGraphGini;
     30    TGraph  fGraphError;
     31    TGraph  fGraphNodes;
     32
     33    TPaveText fRules;
    2534
    2635public:
    2736    MHRanForestGini(Int_t nbins=100, const char *name=NULL, const char *title=NULL);
    28     ~MHRanForestGini();
    29 
    30     TGraph *GetGraphGini() const  { return fGraphGini; }
    3137
    3238    Bool_t SetupFill(const MParList *plist);

  • trunk/MagicSoft/Mars/mranforest/MRanForest.cc

    r7409 r7413  
    166166        if(grid[i]>=grid[i+1])
    167167        {
    168             *fLog<<inf<<"Grid points must be in increasing order! Ignoring grid."<<endl;
     168            *fLog<<warn<<"Grid points must be in increasing order! Ignoring grid."<<endl;
    169169            return;
    170170        }
    … …  
    225225            fClass[j] = int(fHadTrue[j]+0.5);
    226226    }
    227 
    228     return;
    229227}
    230228
    … …  
    272270Bool_t MRanForest::AddTree(MRanTree *rantree=NULL)
    273271{
    274     fRanTree = rantree ? rantree:fRanTree;
     272    fRanTree = rantree ? rantree : fRanTree;
    275273
    276274    if (!fRanTree) return kFALSE;
    … …  
    287285    // access matrix, copy last column (target) preliminarily
    288286    // into fHadTrue
    289     TMatrix mat_tmp = mat->GetM();
    290     int dim         = mat_tmp.GetNcols();
    291     int numdata     = mat_tmp.GetNrows();
    292 
    293     fMatrix=new TMatrix(mat_tmp);
     287    fMatrix = new TMatrix(mat->GetM());
     288
     289    int dim     = fMatrix->GetNcols()-1;
     290    int numdata = fMatrix->GetNrows();
    294291
    295292    fHadTrue.Set(numdata);
    … …  
    297294
    298295    for (Int_t j=0;j<numdata;j++)
    299         fHadTrue[j] = (*fMatrix)(j,dim-1);
     296        fHadTrue[j] = (*fMatrix)(j,dim);
    300297
    301298    // remove last col
    302     fMatrix->ResizeTo(numdata,dim-1);
    303     dim=fMatrix->GetNcols();
     299    fMatrix->ResizeTo(numdata,dim);
    304300
    305301    //-------------------------------------------------------------------
    … …  
    308304    fClass.Reset(0);
    309305
    310     if(fClassify) PrepareClasses();
     306    if (fClassify)
     307        PrepareClasses();
    311308
    312309    //-------------------------------------------------------------------
    313310    // allocating and initializing arrays
    314     fHadEst.Set(numdata);       fHadEst.Reset(0);
    315     fNTimesOutBag.Set(numdata); fNTimesOutBag.Reset(0);
    316     fDataSort.Set(dim*numdata); fDataSort.Reset(0);
    317     fDataRang.Set(dim*numdata); fDataRang.Reset(0);
     311    fHadEst.Set(numdata);
     312    fHadEst.Reset(0);
     313
     314    fNTimesOutBag.Set(numdata);
     315    fNTimesOutBag.Reset(0);
     316
     317    fDataSort.Set(dim*numdata);
     318    fDataSort.Reset(0);
     319
     320    fDataRang.Set(dim*numdata);
     321    fDataRang.Reset(0);
    318322
    319323    if(fWeight.GetSize()!=numdata)
    … …  
    326330    //-------------------------------------------------------------------
    327331    // setup rules to be used for classification/regression
    328     MDataArray *allrules=(MDataArray*)mat->GetColumns();
     332    const MDataArray *allrules=(MDataArray*)mat->GetColumns();
    329333    if(allrules==NULL)
    330334    {
    … …  
    333337    }
    334338
    335     fRules=new MDataArray(); fRules->Reset();
    336     TString target_rule;
    337 
    338     for(Int_t i=0;i<dim+1;i++)
    339     {
    340         MData &data=(*allrules)[i];
    341         if(i<dim)
    342             fRules->AddEntry(data.GetRule());
    343         else
    344             target_rule=data.GetRule();
    345     }
    346 
    347     *fLog << inf <<endl;
    348     *fLog << inf <<"Setting up RF for training on target:"<<endl<<" "<<target_rule.Data()<<endl;
    349     *fLog << inf <<"Following rules are used as input to RF:"<<endl;
    350 
    351     for(Int_t i=0;i<dim;i++)
    352     {
    353         MData &data=(*fRules)[i];
    354         *fLog<<inf<<" "<<i<<") "<<data.GetRule()<<endl<<flush;
    355     }
    356 
    357     *fLog << inf <<endl;
     339    fRules = new MDataArray();
     340    fRules->Reset();
     341
     342    const TString target_rule = (*allrules)[dim];
     343    for (Int_t i=0;i<dim;i++)
     344        fRules->AddEntry((*allrules)[i].GetRule());
     345
     346    *fLog << inf << endl;
     347    *fLog << "Setting up RF for training on target:" << endl;
     348    *fLog << " " << target_rule.Data() << endl;
     349    *fLog << "Following rules are used as input to RF:" << endl;
     350    for (Int_t i=0;i<dim;i++)
     351        *fLog << " " << i << ") " << (*fRules)[i].GetRule() << endl;
     352
     353    *fLog << endl;
    358354
    359355    //-------------------------------------------------------------------
    360356    // prepare (sort) data for fast optimization algorithm
    361     if(!CreateDataSort()) return kFALSE;
     357    if (!CreateDataSort())
     358        return kFALSE;
    362359
    363360    //-------------------------------------------------------------------
    364361    // access and init tree container
    365362    fRanTree = (MRanTree*)plist->FindCreateObj("MRanTree");
    366 
    367363    if(!fRanTree)
    368364    {
    … …  
    371367    }
    372368
     369    const Int_t tryest = TMath::Nint(TMath::Sqrt(dim)+0.5);
     370
     371    *fLog << inf << endl;
     372    *fLog << "Following input for the tree growing are used:"<<endl;
     373    *fLog << " Number of Trees : "<<fNumTrees<<endl;
     374    *fLog << " Number of Trials: "<<(fNumTry==0?tryest:fNumTry)<<(fNumTry==0?" (auto)":"")<<endl;
     375    *fLog << " Final Node size : "<<fNdSize<<endl;
     376    *fLog << " Using Grid:       "<<(fGrid.GetSize()>0?"Yes":"No")<<endl;
     377    *fLog << " Number of Events: "<<numdata<<endl;
     378    *fLog << " Number of Params: "<<dim<<endl;
     379
     380    if(fNumTry==0)
     381    {
     382        fNumTry=tryest;
     383        *fLog << inf << endl;
     384        *fLog << "Set no. of trials to the recommended value of round("<< TMath::Sqrt(dim) <<") = ";
     385        *fLog << fNumTry << endl;
     386
     387
     388    }
     389    fRanTree->SetNumTry(fNumTry);
    373390    fRanTree->SetClassify(fClassify);
    374391    fRanTree->SetNdSize(fNdSize);
    375392
    376     if(fNumTry==0)
    377     {
    378         double ddim = double(dim);
    379 
    380         fNumTry=int(sqrt(ddim)+0.5);
    381         *fLog<<inf<<endl;
    382         *fLog<<inf<<"Set no. of trials to the recommended value of round("<<sqrt(ddim)<<") = ";
    383         *fLog<<inf<<fNumTry<<endl;
    384 
    385     }
    386     fRanTree->SetNumTry(fNumTry);
    387 
    388     *fLog<<inf<<endl;
    389     *fLog<<inf<<"Following settings for the tree growing are used:"<<endl;
    390     *fLog<<inf<<" Number of Trees : "<<fNumTrees<<endl;
    391     *fLog<<inf<<" Number of Trials: "<<fNumTry<<endl;
    392     *fLog<<inf<<" Final Node size : "<<fNdSize<<endl;
    393 
    394393    fTreeNo=0;
    395394
    … …  
    415414    if (fTreeNo==1)
    416415    {
    417         *fLog << inf << endl;
    418         *fLog << underline;
     416        *fLog << inf << endl << underline;
    419417
    420418        if(calcResolution)
    … …  
    435433    TArrayF winbag(numdata); // Initialization includes filling with 0
    436434
    437     float square=0; float mean=0;
     435    float square=0;
     436    float mean=0;
    438437
    439438    for (Int_t n=0; n<numdata; n++)
    … …  
    483482    for (Int_t ievt=0;ievt<numdata;ievt++)
    484483    {
    485         if (jinbag[ievt]>0) continue;
     484        if (jinbag[ievt]>0)
     485            continue;
    486486
    487487        fHadEst[ievt] +=fRanTree->TreeHad((*fMatrix), ievt);
    … …  
    491491
    492492    Int_t n=0;
    493     double ferr=0;
     493    Float_t ferr=0;
    494494
    495495    for (Int_t ievt=0;ievt<numdata;ievt++)
    … …  
    508508    //-------------------------------------------------------------------
    509509    // give running output
    510     *fLog << inf << setw(5)  << fTreeNo;
    511     *fLog << inf << setw(18) << fRanTree->GetNumEndNodes();
    512     *fLog << inf << Form("%18.2f", ferr*100.);
    513     *fLog << inf << endl;
     510    *fLog << setw(5)  << fTreeNo;
     511    *fLog << setw(18) << fRanTree->GetNumEndNodes();
     512    *fLog << Form("%18.2f", ferr*100.);
     513    *fLog << endl;
     514
     515    fRanTree->SetError(ferr);
    514516
    515517    // adding tree to forest

  • trunk/MagicSoft/Mars/mranforest/MRanForestCalc.cc

    r6893 r7413  
    1616!
    1717!
    18 !   Author(s): Thomas Hengstebeck 3/2003 <mailto:hengsteb@alwa02.physik.uni-siegen.de>
    19 !
    20 !   Copyright: MAGIC Software Development, 2000-2003
     18!   Author(s): Thomas Hengstebeck 2/2005 <mailto:hengsteb@physik.hu-berlin.de>
     19!   Author(s): Thomas Bretz 8/2005 <mailto:tbretz@astro.uni-wuerzburg.de>
     20!
     21!   Copyright: MAGIC Software Development, 2000-2005
    2122!
    2223!
    … …  
    2728//  MRanForestCalc
    2829//
    29 //  Calculates the hadroness of an event. It calculates a mean value of all
    30 //  classifications by the trees in a previously grown random forest.
    31 //
    32 //  To use only n trees for your calculation use:
    33 //  MRanForestCalc::SetUseNumTrees(n);
    3430//
    3531////////////////////////////////////////////////////////////////////////////
    … …  
    3834#include <TVector.h>
    3935
    40 #include "MHMatrix.h" // must be before MLogManip.h
    41 #include "MDataArray.h"
     36#include "MHMatrix.h"
    4237
    4338#include "MLog.h"
    4439#include "MLogManip.h"
    4540
     41#include "MData.h"
     42#include "MDataArray.h"
     43
     44#include "MRanForest.h"
     45#include "MParameters.h"
     46
    4647#include "MParList.h"
    47 
    48 #include "MRanTree.h"
    49 #include "MRanForest.h"
    50 
    51 #include "MParameters.h"
    52 
     48#include "MTaskList.h"
    5349#include "MEvtLoop.h"
    54 #include "MTaskList.h"
     50#include "MRanForestGrow.h"
    5551#include "MFillH.h"
    56 #include "MStatusDisplay.h"
    57 #include "MRanForestGrow.h"
    58 #include "MRanForestFill.h"
    59 
    60 #include "MWriteRootFile.h"
    61 #include "MReadTree.h"
    6252
    6353ClassImp(MRanForestCalc);
    … …  
    6555using namespace std;
    6656
    67 static const TString gsDefName  = "MRanForestCalc";
    68 static const TString gsDefTitle = "Tree Classification Loop 1/2";
    69 
    70 // --------------------------------------------------------------------------
    71 //
    72 // Setup histograms and the number of distances which are used for
    73 // avaraging in CalcDist
    74 //
     57const TString MRanForestCalc::gsDefName    = "MRanForestCalc";
     58const TString MRanForestCalc::gsDefTitle   = "RF for energy estimation";
     59
     60const TString MRanForestCalc::gsNameOutput = "RanForestOut";
     61
    7562MRanForestCalc::MRanForestCalc(const char *name, const char *title)
    76     : fNum(100), fHadronnessName("MHadronness"), fData(NULL), fRanForest(0), fRanTree(0)
    77 {
    78     //
    79     //   set the name and title of this object
    80     //
     63    : fDebug(kFALSE), fData(0), fRFOut(0),
     64    fNumTrees(-1), fNumTry(-1), fNdSize(-1), fNumObsoleteVariables(1),
     65    fTestMatrix(0), fEstimationMode(kMean)
     66{
    8167    fName  = name  ? name  : gsDefName.Data();
    8268    fTitle = title ? title : gsDefTitle.Data();
    83 }
    84 
    85 // --------------------------------------------------------------------------
    86 //
    87 // Delete the data chains
    88 //
     69
     70    // FIXME:
     71    fNumTrees = 100; //100
     72    fNumTry   = 0;   //3   0 means: in MRanForest estimated best value will be calculated
     73    fNdSize   = 1;   //1   
     74}
     75
    8976MRanForestCalc::~MRanForestCalc()
    9077{
    91     if (!IsOwner())
    92         return;
    93 
    94     delete fRanForest;
    95     delete fRanTree;
    96 }
    97 
    98 // --------------------------------------------------------------------------
    99 //
    100 // Needs:
    101 //  - MatrixGammas  [MHMatrix]
    102 //  - MatrixHadrons [MHMatrix]
    103 //  - MHadronness
    104 //  - all data containers used to build the matrixes
    105 //
    106 // The matrix object can be filles using MFillH. And must be of the same
    107 // number of columns (with the same meaning).
    108 //
     78    fEForests.Delete();
     79}
     80
     81Int_t MRanForestCalc::Train(const MHMatrix &matrixtrain, const TArrayD &grid, Int_t ver)
     82{
     83    gLog.Separator("MRanForestCalc - Train");
     84
     85    if (!matrixtrain.GetColumns())
     86    {
     87        *fLog << err << "ERROR - MHMatrix does not contain rules... abort." << endl;
     88        return kFALSE;
     89    }
     90
     91    const Int_t ncols = matrixtrain.GetM().GetNcols();
     92    const Int_t nrows = matrixtrain.GetM().GetNrows();
     93    if (ncols<=0 || nrows <=0)
     94    {
     95        *fLog << err << "ERROR - No. of columns or no. of rows of matrixtrain equal 0 ... abort." << endl;
     96        return kFALSE;
     97    }
     98
     99    // rules (= combination of image par) to be used for energy estimation
     100    TFile fileRF(fFileName, "recreate");
     101    if (!fileRF.IsOpen())
     102    {
     103        *fLog << err << "ERROR - File to store RFs could not be opened... abort." << endl;
     104        return kFALSE;
     105    }
     106
     107    const Int_t nobs = fNumObsoleteVariables; // Number of obsolete columns
     108
     109    const MDataArray &dcol = *matrixtrain.GetColumns();
     110
     111    MDataArray usedrules;
     112    for (Int_t i=0; i<ncols; i++)
     113        if (i<ncols-nobs)  // -3 is important!!!
     114            usedrules.AddEntry(dcol[i].GetRule());
     115        else
     116            *fLog << inf << "Skipping " << dcol[i].GetRule() << " for training" << endl;
     117
     118    MDataArray rules(usedrules);
     119    rules.AddEntry(ver<2?"Classification":dcol[ncols-1].GetRule());
     120
     121    // prepare matrix for current energy bin
     122    TMatrix mat(matrixtrain.GetM());
     123
     124    // last column must be removed (true energy col.)
     125    mat.ResizeTo(nrows, ncols-nobs+1);
     126
     127    if (fDebug)
     128        gLog.SetNullOutput(kTRUE);
     129
     130    const Int_t nbins = ver>0 ? 1 : grid.GetSize()-1;
     131    for (Int_t ie=0; ie<nbins; ie++)
     132    {
     133        switch (ver)
     134        {
     135        case 0: // Replace Energy Grid by classification
     136            {
     137                Int_t irows=0;
     138                for (Int_t j=0; j<nrows; j++)
     139                {
     140                    const Double_t energy = matrixtrain.GetM()(j,ncols-1);
     141                    const Bool_t   inside = energy>grid[ie] && energy<=grid[ie+1];
     142
     143                    mat(j, ncols-nobs) = inside ? 1 : 0;
     144
     145                    if (inside)
     146                        irows++;
     147                }
     148                if (irows==0)
     149                    *fLog << warn << "WARNING - Skipping";
     150                else
     151                    *fLog << inf << "Training RF for";
     152
     153                *fLog << " energy bin " << ie << " (" << grid[ie] << ", " << grid[ie+1] << ") " << irows << "/" << nrows << endl;
     154
     155                if (irows==0)
     156                    continue;
     157            }
     158            break;
     159
     160        case 1: // Use Energy as classifier
     161        case 2:
     162            for (Int_t j=0; j<nrows; j++)
     163                mat(j, ncols-nobs) = matrixtrain.GetM()(j,ncols-1);
     164            break;
     165        }
     166
     167        MHMatrix matrix(mat, &rules, "MatrixTrain");
     168
     169        MParList plist;
     170        MTaskList tlist;
     171        plist.AddToList(&tlist);
     172        plist.AddToList(&matrix);
     173
     174        MRanForest rf;
     175        rf.SetNumTrees(fNumTrees);
     176        rf.SetNumTry(fNumTry);
     177        rf.SetNdSize(fNdSize);
     178        rf.SetClassify(ver<2 ? 1 : 0);
     179        if (ver==1)
     180            rf.SetGrid(grid);
     181
     182        plist.AddToList(&rf);
     183
     184        MRanForestGrow rfgrow;
     185        tlist.AddToList(&rfgrow);
     186
     187        MFillH fillh("MHRanForestGini");
     188        tlist.AddToList(&fillh);
     189
     190        MEvtLoop evtloop;
     191        evtloop.SetParList(&plist);
     192        evtloop.SetDisplay(fDisplay);
     193        evtloop.SetLogStream(fLog);
     194
     195        if (!evtloop.Eventloop())
     196            return kFALSE;
     197
     198        if (fDebug)
     199            gLog.SetNullOutput(kFALSE);
     200
     201        if (ver==0)
     202        {
     203            // Calculate bin center
     204            const Double_t E = (TMath::Log10(grid[ie])+TMath::Log10(grid[ie+1]))/2;
     205
     206            // save whole forest
     207            rf.SetUserVal(E);
     208            rf.SetName(Form("%.10f", E));
     209        }
     210
     211        rf.Write();
     212    }
     213
     214    // save rules
     215    usedrules.Write("rules");
     216
     217    return kTRUE;
     218}
     219
     220Int_t MRanForestCalc::ReadForests(MParList &plist)
     221{
     222    TFile fileRF(fFileName, "read");
     223    if (!fileRF.IsOpen())
     224    {
     225        *fLog << err << dbginf << "File containing RFs could not be opened... aborting." << endl;
     226        return kFALSE;
     227    }
     228
     229    fEForests.Delete();
     230
     231    TIter Next(fileRF.GetListOfKeys());
     232    TObject *o=0;
     233    while ((o=Next()))
     234    {
     235        MRanForest *forest=0;
     236        fileRF.GetObject(o->GetName(), forest);
     237        if (!forest)
     238            continue;
     239
     240        forest->SetUserVal(atof(o->GetName()));
     241
     242        fEForests.Add(forest);
     243    }
     244
     245    // Maybe fEForests[0].fRules yould be used instead?
     246
     247    if (fData->Read("rules")<=0)
     248    {
     249        *fLog << err << "ERROR - Reading 'rules' from file " << fFileName << endl;
     250        return kFALSE;
     251    }
     252
     253    return kTRUE;
     254}
     255
    109256Int_t MRanForestCalc::PreProcess(MParList *plist)
    110257{
    111     if (!fRanForest)
    112     {
    113         fRanForest = (MRanForest*)plist->FindObject("MRanForest");
    114         if (!fRanForest)
     258    fRFOut = (MParameterD*)plist->FindCreateObj("MParameterD", fNameOutput);
     259    if (!fRFOut)
     260        return kFALSE;
     261
     262    fData = (MDataArray*)plist->FindCreateObj("MDataArray");
     263    if (!fData)
     264        return kFALSE;
     265
     266    if (!ReadForests(*plist))
     267    {
     268        *fLog << err << "Reading RFs failed... aborting." << endl;
     269        return kFALSE;
     270    }
     271
     272    *fLog << inf << "RF read from " << fFileName << endl;
     273
     274    if (fTestMatrix)
     275        return kTRUE;
     276
     277    fData->Print();
     278
     279    if (!fData->PreProcess(plist))
     280    {
     281        *fLog << err << "PreProcessing of the MDataArray failed... aborting." << endl;
     282        return kFALSE;
     283    }
     284
     285    return kTRUE;
     286}
     287
     288#include <TGraph.h>
     289#include <TF1.h>
     290Int_t MRanForestCalc::Process()
     291{
     292    TVector event;
     293    if (fTestMatrix)
     294        *fTestMatrix >> event;
     295    else
     296        *fData >> event;
     297
     298    // --------------- Single Tree RF -------------------
     299    if (fEForests.GetEntries()==1)
     300    {
     301        MRanForest *rf = (MRanForest*)fEForests[0];
     302        fRFOut->SetVal(rf->CalcHadroness(event));
     303        fRFOut->SetReadyToSave();
     304
     305        return kTRUE;
     306    }
     307
     308    // --------------- Multi Tree RF -------------------
     309    static TF1 f1("f1", "gaus");
     310
     311    Double_t sume = 0;
     312    Double_t sumh = 0;
     313    Double_t maxh = 0;
     314    Double_t maxe = 0;
     315
     316    Double_t max  = -1e10;
     317    Double_t min  =  1e10;
     318
     319    TIter Next(&fEForests);
     320    MRanForest *rf = 0;
     321
     322    TGraph g;
     323    while ((rf=(MRanForest*)Next()))
     324    {
     325        const Double_t h = rf->CalcHadroness(event);
     326        const Double_t e = rf->GetUserVal();
     327
     328        g.SetPoint(g.GetN(), e, h);
     329
     330        sume += e*h;
     331        sumh += h;
     332
     333        if (h>maxh)
    115334        {
    116             *fLog << err << dbginf << "MRanForest not found... aborting." << endl;
    117             return kFALSE;
     335            maxh = h;
     336            maxe = e;
    118337        }
    119     }
    120 
    121     if (!fRanTree)
    122     {
    123         fRanTree = (MRanTree*)plist->FindObject("MRanTree");
    124         if (!fRanTree)
    125         {
    126             *fLog << err << dbginf << "MRanTree not found... aborting." << endl;
    127             return kFALSE;
    128         }
    129     }
    130 
    131     fData = fRanTree->GetRules();
    132 
    133     if (!fData)
    134     {
    135         *fLog << err << dbginf << "Error matrix doesn't contain columns... aborting." << endl;
    136         return kFALSE;
    137     }
    138 
    139     if (!fData->PreProcess(plist))
    140     {
    141         *fLog << err << dbginf << "PreProcessing of the MDataArray failed for the columns failed... aborting." << endl;
    142         return kFALSE;
    143     }
    144 
    145     fHadroness = (MParameterD*)plist->FindCreateObj("MParameterD", fHadronnessName);
    146     if (!fHadroness)
    147         return kFALSE;
     338        if (e>max)
     339            max = e;
     340        if (e<min)
     341            min = e;
     342    }
     343
     344    switch (fEstimationMode)
     345    {
     346    case kMean:
     347        fRFOut->SetVal(pow(10, sume/sumh));
     348        break;
     349    case kMaximum:
     350        fRFOut->SetVal(pow(10, maxe));
     351        break;
     352    case kFit:
     353        f1.SetParameter(0, maxh);
     354        f1.SetParameter(1, maxe);
     355        f1.SetParameter(2, 0.125);
     356        g.Fit(&f1, "Q0N");
     357        fRFOut->SetVal(pow(10, f1.GetParameter(1)));
     358        break;
     359    }
     360
     361    fRFOut->SetReadyToSave();
    148362
    149363    return kTRUE;
    … …  
    153367//
    154368//
    155 Int_t MRanForestCalc::Process()
    156 {
    157     // first copy the data from the data array to a vector event
    158     TVector event;
    159     *fData >> event;
    160 
    161     Double_t hadroness=fRanForest->CalcHadroness(event);
    162     fHadroness->SetVal(hadroness);
    163 
    164     return kTRUE;
    165 }
    166 
    167 Bool_t MRanForestCalc::Grow(MHMatrix *matrixg,MHMatrix *matrixh,Int_t ntree,
    168                             Int_t numtry,Int_t ndsize,const char* treefile,
    169                             const char* treename,const char* contname,
    170                             const char* hgininame)
    171 {
    172 
    173     treename  = treename  ? treename  : "Tree";
    174     contname  = contname  ? contname  : "MRanTree";
    175     hgininame  = hgininame  ? hgininame  : "MHRanForestGini";
    176 
    177     if (!matrixg->IsValid())
    178     {
    179         *fLog << err << dbginf << " MRanForestCalc::Grow - ERROR: matrixg not valid." << endl;
    180         return kFALSE;
    181     }
    182     if(!matrixh->IsValid())
    183     {
    184         *fLog << err << dbginf << " MRanForestCalc::Grow - ERROR: matrixh not valid." << endl;
    185         return kFALSE;
    186     }
    187 
    188     MEvtLoop run(GetName());
    189     MTaskList tlist;
    190     MParList plist;
    191     plist.AddToList(&tlist);
    192     plist.AddToList(matrixg);
    193     plist.AddToList(matrixh);
    194 
    195     // creating training task and setting parameters
    196     MRanForestGrow rfgrow;
    197     rfgrow.SetNumTrees(ntree); // number of trees
    198     rfgrow.SetNumTry(numtry);  // number of trials in random split selection
    199     rfgrow.SetNdSize(ndsize);  // limit for nodesize
    200     tlist.AddToList(&rfgrow);
    201 
    202     if(treefile){
    203         MWriteRootFile rfwrite(treefile);
    204         rfwrite.AddContainer(contname,treename);
    205         tlist.AddToList(&rfwrite);
    206     }
    207 
    208     MFillH fillh(hgininame);
    209     tlist.AddToList(&fillh);
    210 
    211     run.SetParList(&plist);
    212    
    213     // Execute tree growing
    214     if (!run.Eventloop())
    215     {
    216         *fLog << err << dbginf << "Evtloop in MRanForestCalc::Grow failed." << endl;
    217         return kFALSE;
    218     }
    219     tlist.PrintStatistics(0, kTRUE);
    220 
    221     if (TestBit(kEnableGraphicalOutput))
    222         plist.FindObject(hgininame)->DrawClone();
    223 
    224     return kTRUE;
    225 }
    226 
    227 Bool_t MRanForestCalc::Fill(Int_t ntree,const char* treefile,const char* treename)
    228 {
    229     treefile  = treefile  ? treefile  : "RF.root";
    230     treename  = treename  ? treename  : "Tree";
    231 
    232     MParList  plist;
    233 
    234     MTaskList tlist;
    235     plist.AddToList(&tlist);
    236 
    237     MReadTree read(treename,treefile);
    238     read.DisableAutoScheme();
    239 
    240     MRanForestFill rffill;
    241     rffill.SetNumTrees(ntree);
    242 
    243     tlist.AddToList(&read);
    244     tlist.AddToList(&rffill);
    245 
    246     MEvtLoop run(GetName());
    247     run.SetParList(&plist);
    248 
    249     //
    250     // Execute tree reading
    251     //
    252     if (!run.Eventloop())
    253     {
    254         *fLog << err << dbginf << "Evtloop in MRanForestCalc::Fill failed." << endl;
    255         return kFALSE;
    256     }
    257     tlist.PrintStatistics(0, kTRUE);
    258 
    259     return kTRUE;
    260 }
    261 
    262369Int_t MRanForestCalc::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
    263370{
    264     if (!IsEnvDefined(env, prefix, "File",     print))
    265         return kFALSE;
    266 
    267     TString fname = GetEnvValue(env, prefix, "File", "RF.root");
    268     TString tname = GetEnvValue(env, prefix, "Tree", "Tree");
    269     const Int_t   num   = GetEnvValue(env, prefix, "NumTrees", 100);
    270     const Bool_t  debug = GetEnvValue(env, prefix, "Debug", kFALSE);
    271 
    272     fname.ReplaceAll("\015", "");
    273     tname.ReplaceAll("\015", "");
    274 
    275     *fLog << inf << dbginf << "Reading " << num << " trees " << tname << " from file " << fname << endl;
    276 
    277     gLog.SetNullOutput(!debug);
    278     MEvtLoop evtloop;
    279     MParList  plist;
    280     evtloop.SetParList(&plist);
    281     MLog l;
    282     l.SetNullOutput(!debug);
    283     evtloop.SetLogStream(&l);
    284     gLog.SetNullOutput(debug);
    285 
    286     if (IsOwner())
    287     {
    288         delete fRanForest;
    289         delete fRanTree;
    290     }
    291     fRanForest = new MRanForest;
    292     fRanTree   = new MRanTree;
    293     SetOwner();
    294 
    295     plist.AddToList(fRanForest);
    296     plist.AddToList(fRanTree);
    297 
    298     MTaskList tlist;
    299     plist.AddToList(&tlist);
    300 
    301     MReadTree read(tname, fname);
    302     read.DisableAutoScheme();
    303 
    304     MRanForestFill rffill;
    305     rffill.SetNumTrees(num);
    306 
    307     tlist.AddToList(&read);
    308     tlist.AddToList(&rffill);
    309 
    310     if (!evtloop.Eventloop())
    311     {
    312         *fLog << err << "ERROR - Reading " << tname << " from file " << fname << endl;
    313         return kERROR;
    314     }
    315 
    316     return kTRUE;
    317 }
     371    Bool_t rc = kFALSE;
     372    if (IsEnvDefined(env, prefix, "FileName", print))
     373    {
     374        rc = kTRUE;
     375        SetFileName(GetEnvValue(env, prefix, "FileName", fFileName));
     376    }
     377    if (IsEnvDefined(env, prefix, "Debug", print))
     378    {
     379        rc = kTRUE;
     380        SetDebug(GetEnvValue(env, prefix, "Debug", fDebug));
     381    }
     382    if (IsEnvDefined(env, prefix, "NameOutput", print))
     383    {
     384        rc = kTRUE;
     385        SetNameOutput(GetEnvValue(env, prefix, "NameOutput", fNameOutput));
     386    }
     387    if (IsEnvDefined(env, prefix, "EstimationMode", print))
     388    {
     389        TString txt = GetEnvValue(env, prefix, "EstimationMode", "");
     390        txt = txt.Strip(TString::kBoth);
     391        txt.ToLower();
     392        if (txt==(TString)"mean")
     393            fEstimationMode = kMean;
     394        if (txt==(TString)"maximum")
     395            fEstimationMode = kMaximum;
     396        if (txt==(TString)"fit")
     397            fEstimationMode = kFit;
     398        rc = kTRUE;
     399    }
     400    return rc;
     401}

  • trunk/MagicSoft/Mars/mranforest/MRanForestCalc.h

    r6893 r7413  
    66#endif
    77
    8 #ifndef MARS_MHMatrix
    9 #include "MHMatrix.h"
     8#ifndef ROOT_TObjArray
     9#include <TObjArray.h>
    1010#endif
    1111
    12 class MParList;
     12#ifndef ROOT_TArrayD
     13#include <TArrayD.h>
     14#endif
     15
     16class MDataArray;
    1317class MParameterD;
    14 class MDataArray;
    15 class MRanTree;
    16 class MRanForest;
     18class MHMatrix;
    1719
    1820class MRanForestCalc : public MTask
    1921{
     22public:
     23    enum EstimationMode_t
     24    {
     25        kMean,
     26        kMaximum,
     27        kFit
     28    };
     29
    2030private:
    21     Int_t  fNum;              // number of trees used to compute hadronness
     31    static const TString gsDefName;     //! Default Name
     32    static const TString gsDefTitle;    //! Default Title
     33    static const TString gsNameOutput;  //! Default Output name
    2234
    23     TString fHadronnessName;  // Name of container storing hadronness
     35    Bool_t       fDebug;      // Debugging of eventloop while training on/off
    2436
    25     MParameterD *fHadroness;  //! Output container for calculated hadroness
     37    TString      fFileName;   // File name to forest
     38    TObjArray    fEForests;   // List of forests
     39
     40    TString      fNameOutput; // Name of output container
     41
    2642    MDataArray  *fData;       //! Used to store the MDataChains to get the event values
    27     MRanForest  *fRanForest;
    28     MRanTree    *fRanTree;
     43    MParameterD *fRFOut;      //! Used to store result
    2944
     45    Int_t        fNumTrees;   //! Training parameters
     46    Int_t        fNumTry;     //! Training parameters
     47    Int_t        fNdSize;     //! Training parameters
     48
     49    Int_t        fNumObsoleteVariables;
     50
     51    MHMatrix    *fTestMatrix; //! Test Matrix used in Process (together with MMatrixLoop)
     52
     53    EstimationMode_t fEstimationMode;
     54
     55private:
     56    // MTask
    3057    Int_t PreProcess(MParList *plist);
    3158    Int_t Process();
    3259
    33     enum { kIsOwner = BIT(14) };
     60    // MRanForestCalc
     61    Int_t ReadForests(MParList &plist);
    3462
    35     Bool_t IsOwner() const { return TestBit(kIsOwner); }
    36     void SetOwner() { SetBit(kIsOwner); }
     63    // MParContainer
     64    Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print);
     65
     66    // Train Interface
     67    Int_t Train(const MHMatrix &n, const TArrayD &grid, Int_t ver=2);
    3768
    3869public:
    … …  
    4071    ~MRanForestCalc();
    4172
    42     Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print);
     73    // Setter for estimation
     74    void SetFileName(TString filename)            { fFileName = filename; }
     75    void SetEstimationMode(EstimationMode_t op)   { fEstimationMode = op; }
     76    void SetNameOutput(TString name=gsNameOutput) { fNameOutput = name; }
    4377
    44     void SetHadronnessName(const TString name) { fHadronnessName = name; }
    45     TString GetHadronnessName() const { return fHadronnessName; }
     78    // Setter for training
     79    void SetNumTrees(UShort_t n=100) { fNumTrees = n; }
     80    void SetNdSize(UShort_t n=5)     { fNdSize   = n; }
     81    void SetNumTry(UShort_t n=0)     { fNumTry   = n; }
     82    void SetDebug(Bool_t b=kTRUE)    { fDebug    = b; }
    4683
    47     void SetUseNumTrees(UShort_t n=100) { fNum = n; }
     84    void SetNumObsoleteVariables(Int_t n=1) { fNumObsoleteVariables = n; }
    4885
    49     Bool_t Grow(MHMatrix *matrixg,MHMatrix *matrixh,Int_t ntree,
    50                 Int_t numtry,Int_t ndsize,const char* treefile=0,
    51                 const char* treename=0,const char* contname=0,
    52                 const char* hgininame=0);
     86    // Train Interface
     87    Int_t TrainMultiRF(const MHMatrix &n, const TArrayD &grid)
     88    {
     89        return Train(n, grid, 0);
     90    }
     91    Int_t TrainSingleRF(const MHMatrix &n, const TArrayD &grid=TArrayD())
     92    {
     93        return Train(n, grid, grid.GetSize()==0 ? 2 : 1);
     94    }
    5395
    54     Bool_t Fill(Int_t ntree,const char* treefile=0,const char* treename=0);
     96    // Test Interface
     97    void  SetTestMatrix(MHMatrix *m=0) { fTestMatrix=m; }
     98    void  InitMapping(MHMatrix *m=0)   { fTestMatrix=m; }
    5599
    56 
    57     ClassDef(MRanForestCalc, 0) // Task
     100    ClassDef(MRanForestCalc, 0) // Task to calculate RF output and for RF training
    58101};
    59102

  • trunk/MagicSoft/Mars/mranforest/MRanTree.h

    r7396 r7413  
    2828    Int_t fNumNodes;
    2929    Int_t fNumEndNodes;
     30    Float_t fError;
    3031
    3132    TArrayI fBestVar;
    … …  
    7071    void SetNdSize(Int_t n);
    7172    void SetNumTry(Int_t n);
     73    void SetError(Float_t f) { fError = f; }
    7274
    7375    Int_t GetNdSize() const { return fNdSize; }
    … …  
    7577    Int_t GetNumNodes()          const { return fNumNodes; }
    7678    Int_t GetNumEndNodes()       const { return fNumEndNodes; }
     79    Float_t GetError() const { return fError; }
    7780
    7881    Int_t GetBestVar(Int_t i)    const { return fBestVar.At(i); }

  • trunk/MagicSoft/Mars/mranforest/Makefile

    r7396 r7413  
    2525           MRanForest.cc \
    2626           MRanForestGrow.cc \
     27           MRanForestCalc.cc \
    2728           MHRanForest.cc \
    28            MHRanForestGini.cc \
    29            MRFEnergyEst.cc
     29           MHRanForestGini.cc
    3030
    3131############################################################

  • trunk/MagicSoft/Mars/mranforest/RanForestLinkDef.h

    r7396 r7413  
    88#pragma link C++ class MRanForest+;
    99#pragma link C++ class MRanForestGrow+;
     10#pragma link C++ class MRanForestCalc+;
    1011
    1112#pragma link C++ class MHRanForest+;
    1213#pragma link C++ class MHRanForestGini+;
    1314
    14 #pragma link C++ class MRFEnergyEst+;
    15 
    1615#endif
Note: See TracChangeset for help on using the changeset viewer.