tmva » src

PDF.cxx

Go to the documentation of this file.
00001 // @(#)root/tmva $Id: PDF.cxx 37574 2010-12-13 16:50:41Z brun $
00002 // Author: Asen Christov, Andreas Hoecker, Joerg Stelzer, Helge Voss, Kai Voss
00003 
00004 /**********************************************************************************
00005  * Project: TMVA - a Root-integrated toolkit for multivariate data analysis       *
00006  * Package: TMVA                                                                  *
00007  * Class  : PDF                                                                   *
00008  * Web    : http://tmva.sourceforge.net                                           *
00009  *                                                                                *
00010  * Description:                                                                   *
00011  *      Implementation (see header for description)                               *
00012  *                                                                                *
00013  * Authors (alphabetical):                                                        *
00014  *      Asen Christov   <christov@physik.uni-freiburg.de> - Freiburg U., Germany  *
00015  *      Andreas Hoecker <Andreas.Hocker@cern.ch> - CERN, Switzerland              *
00016  *      Helge Voss      <Helge.Voss@cern.ch>     - MPI-K Heidelberg, Germany      *
00017  *      Kai Voss        <Kai.Voss@cern.ch>       - U. of Victoria, Canada         *
00018  *                                                                                *
00019  * Copyright (c) 2005:                                                            *
00020  *      CERN, Switzerland                                                         *
00021  *      U. of Victoria, Canada                                                    *
00022  *      MPI-K Heidelberg, Germany                                                 *
00023  *      Freiburg U., Germany                                                      *
00024  *                                                                                *
00025  * Redistribution and use in source and binary forms, with or without             *
00026  * modification, are permitted according to the terms listed in LICENSE           *
00027  * (http://tmva.sourceforge.net/LICENSE)                                          *
00028  **********************************************************************************/
00029 
00030 #include <iomanip>
00031 #include <cstdlib>
00032 
00033 #include "TMath.h"
00034 #include "TF1.h"
00035 #include "TH1F.h"
00036 #include "TVectorD.h"
00037 #include "Riostream.h"
00038 
00039 #include "TMVA/Tools.h"
00040 #include "TMVA/PDF.h"
00041 #include "TMVA/TSpline1.h"
00042 #include "TMVA/TSpline2.h"
00043 #include "TMVA/Version.h"
00044 
00045 // static configuration settings
00046 const Int_t    TMVA::PDF::fgNbin_PdfHist      = 10000;
00047 const Bool_t   TMVA::PDF::fgManualIntegration = kTRUE;
00048 const Double_t TMVA::PDF::fgEpsilon           = 1.0e-12;
00049 TMVA::PDF*     TMVA::PDF::fgThisPDF           = 0;
00050 
00051 ClassImp(TMVA::PDF)
00052 
00053 //_______________________________________________________________________
00054 TMVA::PDF::PDF( const TString& name, Bool_t norm )
00055    : Configurable   (""),
00056      fUseHistogram  ( kFALSE ),
00057      fPDFName       ( name ),
00058      fNsmooth       ( 0 ),
00059      fMinNsmooth    (-1 ),
00060      fMaxNsmooth    (-1 ),
00061      fNSmoothHist   ( 0 ),
00062      fInterpolMethod( PDF::kSpline2 ),
00063      fSpline        ( 0 ),
00064      fPDFHist       ( 0 ),
00065      fHist          ( 0 ),
00066      fHistOriginal  ( 0 ),
00067      fGraph         ( 0 ),
00068      fIGetVal       ( 0 ),
00069      fHistAvgEvtPerBin  ( 0 ),
00070      fHistDefinedNBins  ( 0 ),
00071      fKDEtypeString     ( 0 ),
00072      fKDEiterString     ( 0 ),
00073      fBorderMethodString( 0 ),
00074      fInterpolateString ( 0 ),
00075      fKDEtype       ( KDEKernel::kNone ),
00076      fKDEiter       ( KDEKernel::kNonadaptiveKDE ),
00077      fKDEborder     ( KDEKernel::kNoTreatment ),
00078      fFineFactor    ( 0. ),
00079      fReadingVersion( 0 ),
00080      fCheckHist     ( kFALSE ),
00081      fNormalize     ( norm ),
00082      fSuffix        ( "" ),
00083      fLogger        ( 0 )
00084 {
00085    // default constructor needed for ROOT I/O
00086    fLogger   = new MsgLogger(this);
00087    fgThisPDF = this;
00088 }
00089 
00090 //_______________________________________________________________________
00091 TMVA::PDF::PDF( const TString& name,
00092                 const TH1 *hist,
00093                 PDF::EInterpolateMethod method,
00094                 Int_t minnsmooth,
00095                 Int_t maxnsmooth,
00096                 Bool_t checkHist,
00097                 Bool_t norm) :
00098    Configurable   (""),
00099    fUseHistogram  ( kFALSE ),
00100    fPDFName       ( name ),
00101    fMinNsmooth    ( minnsmooth ),
00102    fMaxNsmooth    ( maxnsmooth ),
00103    fNSmoothHist   ( 0 ),
00104    fInterpolMethod( method ),
00105    fSpline        ( 0 ),
00106    fPDFHist       ( 0 ),
00107    fHist          ( 0 ),
00108    fHistOriginal  ( 0 ),
00109    fGraph         ( 0 ),
00110    fIGetVal       ( 0 ),
00111    fHistAvgEvtPerBin  ( 0 ),
00112    fHistDefinedNBins  ( 0 ),
00113    fKDEtypeString     ( 0 ),
00114    fKDEiterString     ( 0 ),
00115    fBorderMethodString( 0 ),
00116    fInterpolateString ( 0 ),
00117    fKDEtype       ( KDEKernel::kNone ),
00118    fKDEiter       ( KDEKernel::kNonadaptiveKDE ),
00119    fKDEborder     ( KDEKernel::kNoTreatment ),
00120    fFineFactor    ( 0. ),
00121    fReadingVersion( 0 ),
00122    fCheckHist     ( checkHist ),
00123    fNormalize     ( norm ),
00124    fSuffix        ( "" ),
00125    fLogger        ( 0 )
00126 {
00127    // constructor of spline based PDF:
00128    fLogger   = new MsgLogger(this);
00129    BuildPDF( hist );
00130 }
00131 
00132 //_______________________________________________________________________
00133 TMVA::PDF::PDF( const TString& name,
00134                 const TH1* hist,
00135                 KDEKernel::EKernelType ktype,
00136                 KDEKernel::EKernelIter kiter,
00137                 KDEKernel::EKernelBorder kborder,
00138                 Float_t FineFactor,
00139                 Bool_t norm) :
00140    Configurable   (""),
00141    fUseHistogram  ( kFALSE ),
00142    fPDFName       ( name ),
00143    fNsmooth       ( 0 ),
00144    fMinNsmooth    (-1 ),
00145    fMaxNsmooth    (-1 ),
00146    fNSmoothHist   ( 0 ),
00147    fInterpolMethod( PDF::kKDE ),
00148    fSpline        ( 0 ),
00149    fPDFHist       ( 0 ),
00150    fHist          ( 0 ),
00151    fHistOriginal  ( 0 ),
00152    fGraph         ( 0 ),
00153    fIGetVal       ( 0 ),
00154    fHistAvgEvtPerBin  ( 0 ),
00155    fHistDefinedNBins  ( 0 ),
00156    fKDEtypeString     ( 0 ),
00157    fKDEiterString     ( 0 ),
00158    fBorderMethodString( 0 ),
00159    fInterpolateString ( 0 ),
00160    fKDEtype       ( ktype ),
00161    fKDEiter       ( kiter ),
00162    fKDEborder     ( kborder ),
00163    fFineFactor    ( FineFactor),
00164    fReadingVersion( 0 ),
00165    fCheckHist     ( kFALSE ),
00166    fNormalize     ( norm ),
00167    fSuffix        ( "" ),
00168    fLogger        ( 0 )
00169 {
00170    // constructor of kernel based PDF:
00171    fLogger   = new MsgLogger(this);
00172    BuildPDF( hist );
00173 }
00174 
00175 //_______________________________________________________________________
00176 TMVA::PDF::PDF( const TString& name,
00177                 const TString& options,
00178                 const TString& suffix,
00179                 PDF* defaultPDF,
00180                 Bool_t norm) :
00181    Configurable   (options),
00182    fUseHistogram  ( kFALSE ),
00183    fPDFName       ( name ),
00184    fNsmooth       ( 0 ),
00185    fMinNsmooth    ( -1 ),
00186    fMaxNsmooth    ( -1 ),
00187    fNSmoothHist   ( 0 ),
00188    fInterpolMethod( PDF::kSpline0 ),
00189    fSpline        ( 0 ),
00190    fPDFHist       ( 0 ),
00191    fHist          ( 0 ),
00192    fHistOriginal  ( 0 ),
00193    fGraph         ( 0 ),
00194    fIGetVal       ( 0 ),
00195    fHistAvgEvtPerBin  ( 50 ),
00196    fHistDefinedNBins  ( 0 ),
00197    fKDEtypeString     ( "Gauss" ),
00198    fKDEiterString     ( "Nonadaptive" ),
00199    fBorderMethodString( "None" ),
00200    fInterpolateString ( "Spline2" ),
00201    fKDEtype       ( KDEKernel::kNone ),
00202    fKDEiter       ( KDEKernel::kNonadaptiveKDE ),
00203    fKDEborder     ( KDEKernel::kNoTreatment ),
00204    fFineFactor    ( 1. ),
00205    fReadingVersion( 0 ),
00206    fCheckHist     ( kFALSE ),
00207    fNormalize     ( norm ),
00208    fSuffix        ( suffix ),
00209    fLogger        ( 0 )
00210 {
00211    fLogger   = new MsgLogger(this);
00212    if (defaultPDF != 0) {
00213       fNsmooth            = defaultPDF->fNsmooth;
00214       fMinNsmooth         = defaultPDF->fMinNsmooth;
00215       fMaxNsmooth         = defaultPDF->fMaxNsmooth;
00216       fHistAvgEvtPerBin   = defaultPDF->fHistAvgEvtPerBin;
00217       fHistAvgEvtPerBin   = defaultPDF->fHistAvgEvtPerBin;
00218       fInterpolateString  = defaultPDF->fInterpolateString;
00219       fKDEtypeString      = defaultPDF->fKDEtypeString;
00220       fKDEiterString      = defaultPDF->fKDEiterString;
00221       fFineFactor         = defaultPDF->fFineFactor;
00222       fBorderMethodString = defaultPDF->fBorderMethodString;
00223       fCheckHist          = defaultPDF->fCheckHist;
00224       fHistDefinedNBins   = defaultPDF->fHistDefinedNBins;
00225    }
00226 }
00227 
00228 //___________________fNSmoothHist____________________________________________________
00229 TMVA::PDF::~PDF()
00230 {
00231    // destructor
00232    if (fSpline       != NULL) delete fSpline;
00233    if (fHist         != NULL) delete fHist;
00234    if (fPDFHist      != NULL) delete fPDFHist;
00235    if (fHistOriginal != NULL) delete fHistOriginal;
00236    if (fIGetVal      != NULL) delete fIGetVal;
00237    if (fGraph        != NULL) delete fGraph;
00238    delete fLogger;
00239 }
00240 
00241 //_______________________________________________________________________
00242 void TMVA::PDF::BuildPDF( const TH1* hist )
00243 {
00244    fgThisPDF = this;
00245 
00246    // sanity check
00247    if (hist == NULL) Log() << kFATAL << "Called without valid histogram pointer!" << Endl;
00248 
00249    // histogram should be non empty
00250    if (hist->GetEntries() <= 0)
00251       Log() << kFATAL << "Number of entries <= 0 in histogram: " << hist->GetTitle() << Endl;
00252 
00253    if (fInterpolMethod == PDF::kKDE) {
00254       Log() << "Create "
00255             << ((fKDEiter == KDEKernel::kNonadaptiveKDE) ? "nonadaptive " :
00256                 (fKDEiter == KDEKernel::kAdaptiveKDE)    ? "adaptive "    : "??? ")
00257             << ((fKDEtype == KDEKernel::kGauss)          ? "Gauss "       : "??? ")
00258             << "type KDE kernel for histogram: \"" << hist->GetName() << "\""
00259             << Endl;
00260    }
00261    else {
00262       // another sanity check (nsmooth<0 indicated build with KDE)
00263       if (fMinNsmooth<0)
00264          Log() << kFATAL << "PDF construction called with minnsmooth<0" << Endl;
00265       else if (fMaxNsmooth<=0)
00266          fMaxNsmooth = fMinNsmooth;
00267       else if (fMaxNsmooth<fMinNsmooth)
00268          Log() << kFATAL << "PDF construction called with maxnsmooth<minnsmooth" << Endl;
00269    }
00270 
00271    fHistOriginal = (TH1F*)hist->Clone( TString(hist->GetName()) + "_original" );
00272    fHist         = (TH1F*)hist->Clone( TString(hist->GetName()) + "_smoothed" );
00273    fHistOriginal->SetTitle( fHistOriginal->GetName() ); // reset to new title as well
00274    fHist        ->SetTitle( fHist->GetName() );
00275 
00276    // do not store in current target file
00277    fHistOriginal->SetDirectory(0);
00278    fHist        ->SetDirectory(0);
00279    fUseHistogram = kFALSE;
00280 
00281    if (fInterpolMethod == PDF::kKDE) BuildKDEPDF();
00282    else                              BuildSplinePDF();
00283 }
00284 
00285 //_______________________________________________________________________
00286 Int_t TMVA::PDF::GetHistNBins ( Int_t evtNum )
00287 {
00288    Int_t ResolutionFactor = (fInterpolMethod == PDF::kKDE) ? 5 : 1;
00289    if (evtNum == 0 && fHistDefinedNBins == 0)
00290       Log() << kFATAL << "No number of bins set for PDF" << Endl;
00291    else if (fHistDefinedNBins > 0)
00292       return fHistDefinedNBins * ResolutionFactor;
00293    else if ( evtNum > 0 && fHistAvgEvtPerBin > 0 )
00294       return evtNum / fHistAvgEvtPerBin * ResolutionFactor;
00295    else
00296       Log() << kFATAL << "No number of bins or average event per bin set for PDF" << fHistAvgEvtPerBin << Endl;
00297    return 0;
00298 }
00299 
00300 //_______________________________________________________________________
00301 void TMVA::PDF::BuildSplinePDF()
00302 {
00303    // build the PDF from the original histograms
00304 
00305    // (not useful for discrete distributions, or if no splines are requested)
00306    if (fInterpolMethod != PDF::kSpline0 && fCheckHist) CheckHist();
00307    // use ROOT TH1 smooth method
00308    fNSmoothHist = 0;
00309    if (fMaxNsmooth > 0 && fMinNsmooth >=0 ) SmoothHistogram();
00310 
00311    // fill histogramm to graph
00312    FillHistToGraph();
00313 
00314    //   fGraph->Print();
00315    switch (fInterpolMethod) {
00316 
00317    case kSpline0:
00318       // use original histogram as reference
00319       // this is useful, eg, for discrete variables
00320       fUseHistogram = kTRUE;
00321       break;
00322 
00323    case kSpline1:
00324       fSpline = new TMVA::TSpline1( "spline1", new TGraph(*fGraph) );
00325       break;
00326 
00327    case kSpline2:
00328       fSpline = new TMVA::TSpline2( "spline2", new TGraph(*fGraph) );
00329       break;
00330 
00331    case kSpline3:
00332       fSpline = new TSpline3( "spline3", new TGraph(*fGraph) );
00333       break;
00334 
00335    case kSpline5:
00336       fSpline = new TSpline5( "spline5", new TGraph(*fGraph) );
00337       break;
00338 
00339    default:
00340       Log() << kWARNING << "No valid interpolation method given! Use Spline2" << Endl;
00341       fSpline = new TMVA::TSpline2( "spline2", new TGraph(*fGraph) );
00342       Log() << kFATAL << " Well.. .thinking about it, I better quit so you notice you are forced to fix the mistake " << Endl;
00343       std::exit(1);
00344    }
00345    // fill into histogram
00346    FillSplineToHist();
00347 
00348    if (!UseHistogram()) {
00349       fSpline->SetTitle( (TString)fHist->GetTitle() + fSpline->GetTitle() );
00350       fSpline->SetName ( (TString)fHist->GetName()  + fSpline->GetName()  );
00351    }
00352 
00353 
00354    // sanity check
00355    Double_t integral = GetIntegral();
00356    if (integral < 0) Log() << kFATAL << "Integral: " << integral << " <= 0" << Endl;
00357 
00358    // normalize
00359    if (fNormalize)
00360       if (integral>0) fPDFHist->Scale( 1.0/integral );
00361 
00362    fPDFHist->SetDirectory(0);
00363 }
00364 
00365 //_______________________________________________________________________
00366 void TMVA::PDF::BuildKDEPDF()
00367 {
00368    // creates high-binned reference histogram to be used instead of the
00369    // PDF for speed reasons
00370    fPDFHist = new TH1F( "", "", fgNbin_PdfHist, GetXmin(), GetXmax() );
00371    fPDFHist->SetTitle( (TString)fHist->GetTitle() + "_hist from_KDE" );
00372    fPDFHist->SetName ( (TString)fHist->GetName()  + "_hist_from_KDE" );
00373 
00374    // create the kernel object
00375    Float_t histoLowEdge   = fHist->GetBinLowEdge(1);
00376    Float_t histoUpperEdge = fPDFHist->GetBinLowEdge(fPDFHist->GetNbinsX()) + fPDFHist->GetBinWidth(fPDFHist->GetNbinsX());
00377    KDEKernel *kern = new TMVA::KDEKernel( fKDEiter,
00378                                           fHist, histoLowEdge, histoUpperEdge,
00379                                           fKDEborder, fFineFactor );
00380    kern->SetKernelType(fKDEtype);
00381 
00382    for (Int_t i=1;i<fHist->GetNbinsX();i++) {
00383       // loop over the bins of the original histo
00384       for (Int_t j=1;j<fPDFHist->GetNbinsX();j++) {
00385          // loop over the bins of the new histo and fill it
00386          fPDFHist->AddBinContent(j,fHist->GetBinContent(i)*
00387                                  kern->GetBinKernelIntegral(fPDFHist->GetBinLowEdge(j),
00388                                                             fPDFHist->GetBinLowEdge(j+1),
00389                                                             fHist->GetBinCenter(i),
00390                                                             i)
00391                                  );
00392       }
00393       if (fKDEborder == 3) { // mirror the saples and fill them again
00394          // in order to save time do the mirroring only for the first (the lowwer) 1/5 of the histo to the left;
00395          // and the last (the higher) 1/5 of the histo to the right.
00396          // the middle part of the histo, which is not mirrored, has no influence on the border effects anyway ...
00397          if (i < fHist->GetNbinsX()/5  ) {  // the first (the lowwer) 1/5 of the histo
00398             for (Int_t j=1;j<fPDFHist->GetNbinsX();j++) {
00399                // loop over the bins of the PDF histo and fill it
00400                fPDFHist->AddBinContent(j,fHist->GetBinContent(i)*
00401                                        kern->GetBinKernelIntegral(fPDFHist->GetBinLowEdge(j),
00402                                                                   fPDFHist->GetBinLowEdge(j+1),
00403                                                                   2*histoLowEdge-fHist->GetBinCenter(i), //  mirroring to the left
00404                                                                   i)
00405                                        );
00406             }
00407          }
00408          if (i > 4*fHist->GetNbinsX()/5) { // the last (the higher) 1/5 of the histo
00409             for (Int_t j=1;j<fPDFHist->GetNbinsX();j++) {
00410                // loop over the bins of the PDF histo and fill it
00411                fPDFHist->AddBinContent( j,fHist->GetBinContent(i)*
00412                                         kern->GetBinKernelIntegral(fPDFHist->GetBinLowEdge(j),
00413                                                                    fPDFHist->GetBinLowEdge(j+1),
00414                                                                    2*histoUpperEdge-fHist->GetBinCenter(i), i) );
00415             }
00416          }
00417       }
00418    }
00419 
00420    fPDFHist->SetEntries(fHist->GetEntries());
00421 
00422    delete kern;
00423 
00424    // sanity check
00425    Double_t integral = GetIntegral();
00426    if (integral < 0) Log() << kFATAL << "Integral: " << integral << " <= 0" << Endl;
00427 
00428    // normalize
00429    if (fNormalize)
00430       if (integral>0) fPDFHist->Scale( 1.0/integral );
00431    fPDFHist->SetDirectory(0);
00432 }
00433 
00434 //_______________________________________________________________________
00435 void TMVA::PDF::SmoothHistogram()
00436 {
00437    if(fHist->GetNbinsX()==1) return;
00438    if (fMaxNsmooth == fMinNsmooth) {
00439       fHist->Smooth( fMinNsmooth );
00440       return;
00441    }
00442 
00443    //calculating Mean, RMS of the relative errors and using them to set
00444    //the bounderies of the liniar transformation
00445    Float_t Err=0, ErrAvg=0, ErrRMS=0 ; Int_t num=0, smooth;
00446    for (Int_t bin=0; bin<fHist->GetNbinsX(); bin++) {
00447       if (fHist->GetBinContent(bin+1) <= fHist->GetBinError(bin+1)) continue;
00448       Err = fHist->GetBinError(bin+1) / fHist->GetBinContent(bin+1);
00449       ErrAvg += Err; ErrRMS += Err*Err; num++;
00450    }
00451    ErrAvg /= num;
00452    ErrRMS = TMath::Sqrt(ErrRMS/num-ErrAvg*ErrAvg) ;
00453 
00454    //liniarly convent the histogram to a vector of smoothnum
00455    Float_t MaxErr=ErrAvg+ErrRMS, MinErr=ErrAvg-ErrRMS;
00456    fNSmoothHist = new TH1I("","",fHist->GetNbinsX(),0,fHist->GetNbinsX());
00457    fNSmoothHist->SetTitle( (TString)fHist->GetTitle() + "_Nsmooth" );
00458    fNSmoothHist->SetName ( (TString)fHist->GetName()  + "_Nsmooth" );
00459    for (Int_t bin=0; bin<fHist->GetNbinsX(); bin++) {
00460       if (fHist->GetBinContent(bin+1) <= fHist->GetBinError(bin+1))
00461          smooth=fMaxNsmooth;
00462       else {
00463          Err = fHist->GetBinError(bin+1) / fHist->GetBinContent(bin+1);
00464          smooth=(Int_t)((Err-MinErr) /(MaxErr-MinErr) * (fMaxNsmooth-fMinNsmooth)) + fMinNsmooth;
00465       }
00466       smooth = TMath::Max(fMinNsmooth,TMath::Min(fMaxNsmooth,smooth));
00467       fNSmoothHist->SetBinContent(bin+1,smooth);
00468    }
00469 
00470    //find regions of constant smoothnum, starting from the highest amount of
00471    //smoothing. So the last iteration all the histogram will be smoothed as a whole
00472    for (Int_t n=fMaxNsmooth; n>=0; n--) {
00473       //all the histogram has to be smoothed at least fMinNsmooth
00474       if (n <= fMinNsmooth) { fHist->Smooth(); continue; }
00475       Int_t MinBin=-1,MaxBin =-1;
00476       for (Int_t bin=0; bin < fHist->GetNbinsX(); bin++) {
00477          if (fNSmoothHist->GetBinContent(bin+1) >= n) {
00478             if (MinBin==-1) MinBin = bin;
00479             else MaxBin=bin;
00480          }
00481          else if (MaxBin >= 0) {
00482 #if ROOT_VERSION_CODE > ROOT_VERSION(5,19,2)
00483             fHist->Smooth(1,"R");
00484 #else
00485             fHist->Smooth(1,MinBin+1,MaxBin+1);
00486 #endif
00487             MinBin=MaxBin=-1;
00488          }
00489          else     // can't smooth a single bin
00490             MinBin = -1;
00491       }
00492    }
00493 }
00494 
00495 //_______________________________________________________________________
00496 void TMVA::PDF::FillHistToGraph()
00497 {
00498    // Simple conversion
00499    fGraph=new TGraph(fHist);
00500 }
00501 
00502 //_______________________________________________________________________
00503 void TMVA::PDF::FillSplineToHist()
00504 {
00505    // creates high-binned reference histogram to be used instead of the
00506    // PDF for speed reasons
00507 
00508    if (UseHistogram()) {
00509       // no spline given, use the original histogram
00510       fPDFHist = (TH1*)fHist->Clone();
00511       fPDFHist->SetTitle( (TString)fHist->GetTitle() + "_hist from_spline0" );
00512       fPDFHist->SetName ( (TString)fHist->GetName()  + "_hist_from_spline0" );
00513    }
00514    else {
00515       // create new reference histogram
00516       fPDFHist = new TH1F( "", "", fgNbin_PdfHist, GetXmin(), GetXmax() );
00517       fPDFHist->SetTitle( (TString)fHist->GetTitle() + "_hist from_" + fSpline->GetTitle() );
00518       fPDFHist->SetName ( (TString)fHist->GetName()  + "_hist_from_" + fSpline->GetTitle() );
00519 
00520       for (Int_t bin=1; bin <= fgNbin_PdfHist; bin++) {
00521          Double_t x = fPDFHist->GetBinCenter( bin );
00522          Double_t y = fSpline->Eval( x );
00523          // sanity correction: in cases where strong slopes exist, accidentally, the
00524          // splines can go to zero; in this case we set the corresponding bin content
00525          // equal to the bin content of the original histogram
00526          if (y <= fgEpsilon) y = fHist->GetBinContent( fHist->FindBin( x ) );
00527          fPDFHist->SetBinContent( bin, TMath::Max(y, fgEpsilon) );
00528       }
00529    }
00530    fPDFHist->SetDirectory(0);
00531 }
00532 
00533 //_______________________________________________________________________
00534 void TMVA::PDF::CheckHist() const
00535 {
00536    // sanity check: compare PDF with original histogram
00537    if (fHist == NULL) {
00538       Log() << kFATAL << "<CheckHist> Called without valid histogram pointer!" << Endl;
00539    }
00540 
00541    Int_t nbins = fHist->GetNbinsX();
00542 
00543    Int_t emptyBins=0;
00544    // count number of empty bins
00545    for (Int_t bin=1; bin<=nbins; bin++)
00546       if (fHist->GetBinContent(bin) == 0) emptyBins++;
00547 
00548    if (((Float_t)emptyBins/(Float_t)nbins) > 0.5) {
00549       Log() << kWARNING << "More than 50% (" << (((Float_t)emptyBins/(Float_t)nbins)*100)
00550             <<"%) of the bins in hist '"
00551             << fHist->GetName() << "' are empty!" << Endl;
00552       Log() << kWARNING << "X_min=" << GetXmin()
00553             << " mean=" << fHist->GetMean() << " X_max= " << GetXmax() << Endl;
00554    }
00555 }
00556 
00557 //_______________________________________________________________________
00558 void TMVA::PDF::ValidatePDF( TH1* originalHist ) const
00559 {
00560    // comparison of original histogram with reference PDF
00561 
00562    // if no histogram is given, use the original one (the one the PDF was made from)
00563    if (!originalHist) originalHist = fHistOriginal;
00564 
00565    Int_t    nbins = originalHist->GetNbinsX();
00566 
00567    // treat errors properly
00568    if (originalHist->GetSumw2()->GetSize() == 0) originalHist->Sumw2();
00569 
00570    // ---- first validation: simple(st) possible chi2 test
00571    // count number of empty bins
00572    Double_t chi2 = 0;
00573    Int_t    ndof = 0;
00574    Int_t    nc1  = 0; // deviation counters
00575    Int_t    nc2  = 0; // deviation counters
00576    Int_t    nc3  = 0; // deviation counters
00577    Int_t    nc6  = 0; // deviation counters
00578    for (Int_t bin=1; bin<=nbins; bin++) {
00579       Double_t x  = originalHist->GetBinCenter( bin );
00580       Double_t y  = originalHist->GetBinContent( bin );
00581       Double_t ey = originalHist->GetBinError( bin );
00582 
00583       Int_t binPdfHist = fPDFHist->FindBin( x );
00584       if (binPdfHist<0) continue; // happens only if hist-dim>3
00585 
00586       Double_t yref = GetVal( x );
00587       Double_t rref = ( originalHist->GetSumOfWeights()/fPDFHist->GetSumOfWeights() *
00588                         originalHist->GetBinWidth( bin )/fPDFHist->GetBinWidth( binPdfHist ) );
00589 
00590       if (y > 0) {
00591          ndof++;
00592          Double_t d = TMath::Abs( (y - yref*rref)/ey );
00593          //         cout << "bin: " << bin << "  val: " << x << "  data(err): " << y << "(" << ey << ")   pdf: " 
00594          //              << yref << "  dev(chi2): " << d << "(" << chi2 << ")  rref: " << rref << endl;
00595          chi2 += d*d;
00596          if (d > 1) { nc1++; if (d > 2) { nc2++; if (d > 3) { nc3++; if (d > 6) nc6++; } } }
00597       }
00598    }
00599 
00600    Log() << "Validation result for PDF \"" << originalHist->GetTitle() << "\"" << ": " << Endl;
00601    Log() << Form( "    chi2/ndof(!=0) = %.1f/%i = %.2f (Prob = %.2f)",
00602                   chi2, ndof, chi2/ndof, TMath::Prob( chi2, ndof ) ) << Endl;
00603    if ((1.0 - TMath::Prob( chi2, ndof )) > 0.9999994) {
00604       Log() << kWARNING << "Comparison of the original histogram \"" << originalHist->GetTitle() << "\"" << Endl;
00605       Log() << kWARNING << "with the corresponding PDF gave a chi2/ndof of " << chi2/ndof << "," << Endl;
00606       Log() << kWARNING << "which corresponds to a deviation of more than 5 sigma! Please check!" << Endl;
00607    }
00608    Log() << Form( "    #bins-found(#expected-bins) deviating > [1,2,3,6] sigmas: " \
00609                   "[%i(%i),%i(%i),%i(%i),%i(%i)]",
00610                   nc1, Int_t(TMath::Prob(1.0,1)*ndof), nc2, Int_t(TMath::Prob(4.0,1)*ndof),
00611                   nc3, Int_t(TMath::Prob(9.0,1)*ndof), nc6, Int_t(TMath::Prob(36.0,1)*ndof) ) << Endl;
00612 }
00613 
00614 //_______________________________________________________________________
00615 Double_t TMVA::PDF::GetIntegral() const
00616 {
00617    // computes normalisation
00618 
00619    Double_t integral = fPDFHist->GetSumOfWeights();
00620    integral *= GetPdfHistBinWidth();
00621 
00622    return integral;
00623 }
00624 
00625 //_______________________________________________________________________
00626 Double_t TMVA::PDF::IGetVal( Double_t* x, Double_t* )
00627 {
00628    // static external auxiliary function (integrand)
00629    return ThisPDF()->GetVal( x[0] );
00630 }
00631 
00632 //_______________________________________________________________________
00633 Double_t TMVA::PDF::GetIntegral( Double_t xmin, Double_t xmax )
00634 {
00635    // computes PDF integral within given ranges
00636    Double_t  integral = 0;
00637 
00638    if (fgManualIntegration) {
00639 
00640       // compute integral by summing over bins
00641       Int_t imin = fPDFHist->FindBin(xmin);
00642       Int_t imax = fPDFHist->FindBin(xmax);
00643       if (imin < 1)                     imin = 1;
00644       if (imax > fPDFHist->GetNbinsX()) imax = fPDFHist->GetNbinsX();
00645 
00646       for (Int_t bini = imin; bini <= imax; bini++) {
00647          Float_t dx = fPDFHist->GetBinWidth(bini);
00648          // correct for bin fractions in first and last bin
00649          if      (bini == imin) dx = fPDFHist->GetBinLowEdge(bini+1) - xmin;
00650          else if (bini == imax) dx = xmax - fPDFHist->GetBinLowEdge(bini);
00651          if (dx < 0 && dx > -1.0e-8) dx = 0; // protect against float->double numerical effects
00652          if (dx<0) {
00653             Log() << kWARNING
00654                   << "dx   = " << dx   << std::endl
00655                   << "bini = " << bini << std::endl
00656                   << "xmin = " << xmin << std::endl
00657                   << "xmax = " << xmax << std::endl
00658                   << "imin = " << imin << std::endl
00659                   << "imax = " << imax << std::endl
00660                   << "low edge of imin" << fPDFHist->GetBinLowEdge(imin) << std::endl
00661                   << "low edge of imin+1" << fPDFHist->GetBinLowEdge(imin+1) << Endl;
00662             Log() << kFATAL << "<GetIntegral> dx = " << dx << " < 0" << Endl;
00663          }
00664          integral += fPDFHist->GetBinContent(bini)*dx;
00665       }
00666 
00667    }
00668    else {
00669 
00670       // compute via Gaussian quadrature (C++ version of CERNLIB function DGAUSS)
00671       if (fIGetVal == 0) fIGetVal = new TF1( "IGetVal", PDF::IGetVal, GetXmin(), GetXmax(), 0 );
00672       integral = fIGetVal->Integral( xmin, xmax );
00673    }
00674 
00675    return integral;
00676 }
00677 
00678 //_____________________________________________________________________
00679 Double_t TMVA::PDF::GetVal( Double_t x ) const
00680 {
00681    // returns value PDF(x)
00682 
00683    // check which is filled
00684    Int_t bin = fPDFHist->FindBin(x);
00685    bin = TMath::Max(bin,1);
00686    bin = TMath::Min(bin,fPDFHist->GetNbinsX());
00687 
00688    Double_t retval = 0;
00689 
00690    if (UseHistogram()) {
00691       // use directly histogram bins (this is for discrete PDFs)
00692       retval = fPDFHist->GetBinContent( bin );
00693    }
00694    else {
00695       // linear interpolation
00696       Int_t nextbin = bin;
00697       if ((x > fPDFHist->GetBinCenter(bin) && bin != fPDFHist->GetNbinsX()) || bin == 1)
00698          nextbin++;
00699       else
00700          nextbin--;
00701 
00702       // linear interpolation between adjacent bins
00703       Double_t dx = fPDFHist->GetBinCenter( bin )  - fPDFHist->GetBinCenter( nextbin );
00704       Double_t dy = fPDFHist->GetBinContent( bin ) - fPDFHist->GetBinContent( nextbin );
00705       retval = fPDFHist->GetBinContent( bin ) + (x - fPDFHist->GetBinCenter( bin ))*dy/dx;
00706    }
00707 
00708    return TMath::Max( retval, fgEpsilon );
00709 }
00710 
00711 //_______________________________________________________________________
00712 void TMVA::PDF::DeclareOptions()
00713 {
00714    // define the options (their key words) that can be set in the option string
00715    // know options:
00716    // PDFInterpol[ivar] <string>   Spline0, Spline1, Spline2 <default>, Spline3, Spline5, KDE  used to interpolate reference histograms
00717    //             if no variable index is given, it is valid for ALL the variables
00718    //
00719    // NSmooth           <int>    how often the input histos are smoothed
00720    // MinNSmooth        <int>    min number of smoothing iterations, for bins with most data
00721    // MaxNSmooth        <int>    max number of smoothing iterations, for bins with least data
00722    // NAvEvtPerBin      <int>    minimum average number of events per PDF bin
00723    // TransformOutput   <bool>   transform (often strongly peaked) likelihood output through sigmoid inversion
00724    // fKDEtype          <KernelType>   type of the Kernel to use (1 is Gaussian)
00725    // fKDEiter          <KerneIter>    number of iterations (1 --> "static KDE", 2 --> "adaptive KDE")
00726    // fBorderMethod     <KernelBorder> the method to take care about "border" effects (1=no treatment , 2=kernel renormalization, 3=sample mirroring)
00727 
00728    DeclareOptionRef( fNsmooth, Form("NSmooth%s",fSuffix.Data()),
00729                      "Number of smoothing iterations for the input histograms" );
00730    DeclareOptionRef( fMinNsmooth, Form("MinNSmooth%s",fSuffix.Data()),
00731                      "Min number of smoothing iterations, for bins with most data" );
00732 
00733    DeclareOptionRef( fMaxNsmooth, Form("MaxNSmooth%s",fSuffix.Data()),
00734                      "Max number of smoothing iterations, for bins with least data" );
00735 
00736    DeclareOptionRef( fHistAvgEvtPerBin, Form("NAvEvtPerBin%s",fSuffix.Data()),
00737                      "Average number of events per PDF bin" );
00738 
00739    DeclareOptionRef( fHistDefinedNBins, Form("Nbins%s",fSuffix.Data()),
00740                      "Defined number of bins for the histogram from which the PDF is created" );
00741 
00742    DeclareOptionRef( fCheckHist, Form("CheckHist%s",fSuffix.Data()),
00743                      "Whether or not to check the source histogram of the PDF" );
00744 
00745    DeclareOptionRef( fInterpolateString, Form("PDFInterpol%s",fSuffix.Data()),
00746                      "Interpolation method for reference histograms (e.g. Spline2 or KDE)" );
00747    AddPreDefVal(TString("Spline0")); // take histogram
00748    AddPreDefVal(TString("Spline1")); // linear interpolation between bins
00749    AddPreDefVal(TString("Spline2")); // quadratic interpolation
00750    AddPreDefVal(TString("Spline3")); // cubic interpolation
00751    AddPreDefVal(TString("Spline5")); // fifth order polynome interpolation
00752    AddPreDefVal(TString("KDE"));     // use kernel density estimator
00753 
00754    DeclareOptionRef( fKDEtypeString, Form("KDEtype%s",fSuffix.Data()), "KDE kernel type (1=Gauss)" );
00755    AddPreDefVal(TString("Gauss"));
00756 
00757    DeclareOptionRef( fKDEiterString, Form("KDEiter%s",fSuffix.Data()), "Number of iterations (1=non-adaptive, 2=adaptive)" );
00758    AddPreDefVal(TString("Nonadaptive"));
00759    AddPreDefVal(TString("Adaptive"));
00760 
00761    DeclareOptionRef( fFineFactor , Form("KDEFineFactor%s",fSuffix.Data()),
00762                      "Fine tuning factor for Adaptive KDE: Factor to multyply the width of the kernel");
00763 
00764    DeclareOptionRef( fBorderMethodString, Form("KDEborder%s",fSuffix.Data()),
00765                      "Border effects treatment (1=no treatment , 2=kernel renormalization, 3=sample mirroring)" );
00766    AddPreDefVal(TString("None"));
00767    AddPreDefVal(TString("Renorm"));
00768    AddPreDefVal(TString("Mirror"));
00769 
00770    SetConfigName( GetName() );
00771    SetConfigDescription( "Configuration options for the PDF class" );
00772 }
00773 
00774 //_______________________________________________________________________
00775 void TMVA::PDF::ProcessOptions()
00776 {
00777 
00778    if (fNsmooth < 0) fNsmooth = 0; // set back to default
00779 
00780    if (fMaxNsmooth < 0 || fMinNsmooth < 0) { // use "Nsmooth" variable
00781       fMinNsmooth = fMaxNsmooth = fNsmooth;
00782    }
00783 
00784    if (fMaxNsmooth < fMinNsmooth && fMinNsmooth >= 0) { // sanity check
00785       Log() << kFATAL << "ERROR: MaxNsmooth = "
00786             << fMaxNsmooth << " < MinNsmooth = " << fMinNsmooth << Endl;
00787    }
00788 
00789    if (fMaxNsmooth < 0 || fMinNsmooth < 0) {
00790       Log() << kFATAL << "ERROR: MaxNsmooth = "
00791             << fMaxNsmooth << " or MinNsmooth = " << fMinNsmooth << " smaller than zero" << Endl;
00792    }
00793 
00794    //processing the options
00795    if      (fInterpolateString == "Spline0") fInterpolMethod = TMVA::PDF::kSpline0;
00796    else if (fInterpolateString == "Spline1") fInterpolMethod = TMVA::PDF::kSpline1;
00797    else if (fInterpolateString == "Spline2") fInterpolMethod = TMVA::PDF::kSpline2;
00798    else if (fInterpolateString == "Spline3") fInterpolMethod = TMVA::PDF::kSpline3;
00799    else if (fInterpolateString == "Spline5") fInterpolMethod = TMVA::PDF::kSpline5;
00800    else if (fInterpolateString == "KDE"    ) fInterpolMethod = TMVA::PDF::kKDE;
00801    else if (fInterpolateString != ""       ) {
00802       Log() << kFATAL << "unknown setting for option 'InterpolateMethod': " << fKDEtypeString << ((fSuffix=="")?"":Form(" for pdf with suffix %s",fSuffix.Data())) << Endl;
00803    }
00804 
00805    // init KDE options
00806    if      (fKDEtypeString == "Gauss"      ) fKDEtype = KDEKernel::kGauss;
00807    else if (fKDEtypeString != ""           )
00808       Log() << kFATAL << "unknown setting for option 'KDEtype': " << fKDEtypeString << ((fSuffix=="")?"":Form(" for pdf with suffix %s",fSuffix.Data())) << Endl;
00809    if      (fKDEiterString == "Nonadaptive") fKDEiter = KDEKernel::kNonadaptiveKDE;
00810    else if (fKDEiterString == "Adaptive"   ) fKDEiter = KDEKernel::kAdaptiveKDE;
00811    else if (fKDEiterString != ""           )// nothing more known
00812       Log() << kFATAL << "unknown setting for option 'KDEiter': " << fKDEtypeString << ((fSuffix=="")?"":Form(" for pdf with suffix %s",fSuffix.Data())) << Endl;
00813 
00814    if       ( fBorderMethodString == "None"   ) fKDEborder= KDEKernel::kNoTreatment;
00815    else if  ( fBorderMethodString == "Renorm" ) fKDEborder= KDEKernel::kKernelRenorm;
00816    else if  ( fBorderMethodString == "Mirror" ) fKDEborder= KDEKernel::kSampleMirror;
00817    else if  ( fKDEiterString != ""            ) { // nothing more known
00818       Log() << kFATAL << "unknown setting for option 'KDEBorder': " << fKDEtypeString << ((fSuffix=="")?"":Form(" for pdf with suffix %s",fSuffix.Data())) << Endl;
00819    }
00820 }
00821 
00822 //_______________________________________________________________________
00823 void TMVA::PDF::AddXMLTo( void* parent )
00824 {
00825    // XML file writing
00826    void* pdfxml = gTools().AddChild(parent, "PDF");
00827    gTools().AddAttr(pdfxml, "Name",           fPDFName );
00828    gTools().AddAttr(pdfxml, "MinNSmooth",     fMinNsmooth );
00829    gTools().AddAttr(pdfxml, "MaxNSmooth",     fMaxNsmooth );
00830    gTools().AddAttr(pdfxml, "InterpolMethod", fInterpolMethod );
00831    gTools().AddAttr(pdfxml, "KDE_type",       fKDEtype );
00832    gTools().AddAttr(pdfxml, "KDE_iter",       fKDEiter );
00833    gTools().AddAttr(pdfxml, "KDE_border",     fKDEborder );
00834    gTools().AddAttr(pdfxml, "KDE_finefactor", fFineFactor );
00835    void* pdfhist = gTools().AddChild(pdfxml,"Histogram" );
00836    TH1*  histToWrite = GetOriginalHist();
00837    Bool_t hasEquidistantBinning = gTools().HistoHasEquidistantBins(*histToWrite);
00838    gTools().AddAttr(pdfhist, "Name",  histToWrite->GetName() );
00839    gTools().AddAttr(pdfhist, "NBins", histToWrite->GetNbinsX() );
00840    gTools().AddAttr(pdfhist, "XMin",  histToWrite->GetXaxis()->GetXmin() );
00841    gTools().AddAttr(pdfhist, "XMax",  histToWrite->GetXaxis()->GetXmax() );
00842    gTools().AddAttr(pdfhist, "HasEquidistantBins", hasEquidistantBinning );
00843 
00844    TString bincontent("");
00845    for (Int_t i=0; i<histToWrite->GetNbinsX(); i++) {
00846       bincontent += gTools().StringFromDouble(histToWrite->GetBinContent(i+1));
00847       bincontent += " ";
00848    }
00849    gTools().AddRawLine(pdfhist, bincontent );
00850 
00851    if (!hasEquidistantBinning) {
00852       void* pdfhistbins = gTools().AddChild(pdfxml,"HistogramBinning" );
00853       gTools().AddAttr(pdfhistbins, "NBins", histToWrite->GetNbinsX() );
00854       TString binns("");
00855       for (Int_t i=1; i<=histToWrite->GetNbinsX()+1; i++) {
00856          binns += gTools().StringFromDouble(histToWrite->GetXaxis()->GetBinLowEdge(i));
00857          binns += " ";
00858       }
00859       gTools().AddRawLine(pdfhistbins, binns );
00860    }
00861 }
00862 
00863 //_______________________________________________________________________
00864 void TMVA::PDF::ReadXML( void* pdfnode )
00865 {
00866    // XML file reading
00867    UInt_t enumint;
00868 
00869    gTools().ReadAttr(pdfnode, "MinNSmooth",     fMinNsmooth );
00870    gTools().ReadAttr(pdfnode, "MaxNSmooth",     fMaxNsmooth );
00871    gTools().ReadAttr(pdfnode, "InterpolMethod", enumint ); fInterpolMethod = EInterpolateMethod(enumint);
00872    gTools().ReadAttr(pdfnode, "KDE_type",       enumint ); fKDEtype        = KDEKernel::EKernelType(enumint);
00873    gTools().ReadAttr(pdfnode, "KDE_iter",       enumint ); fKDEiter        = KDEKernel::EKernelIter(enumint);
00874    gTools().ReadAttr(pdfnode, "KDE_border",     enumint ); fKDEborder      = KDEKernel::EKernelBorder(enumint);
00875    gTools().ReadAttr(pdfnode, "KDE_finefactor", fFineFactor );
00876 
00877    TString hname;
00878    UInt_t nbins;
00879    Double_t xmin, xmax;
00880    Bool_t hasEquidistantBinning;
00881 
00882    void* histch = gTools().GetChild(pdfnode);
00883    gTools().ReadAttr( histch, "Name",  hname );
00884    gTools().ReadAttr( histch, "NBins", nbins );
00885    gTools().ReadAttr( histch, "XMin",  xmin );
00886    gTools().ReadAttr( histch, "XMax",  xmax );
00887    gTools().ReadAttr( histch, "HasEquidistantBins", hasEquidistantBinning );
00888 
00889    // recreate the original hist
00890    TH1* newhist = 0;
00891    if (hasEquidistantBinning) {
00892       newhist = new TH1F( hname, hname, nbins, xmin, xmax );
00893       newhist->SetDirectory(0);
00894       const char* content = gTools().GetContent(histch);
00895       std::stringstream s(content);
00896       Double_t val;
00897       for (UInt_t i=0; i<nbins; i++) {
00898          s >> val;
00899          newhist->SetBinContent(i+1,val);
00900       }
00901    }
00902    else {
00903       const char* content = gTools().GetContent(histch);
00904       std::stringstream s(content);
00905       Double_t val;
00906       void* binch = gTools().GetNextChild(histch);
00907       UInt_t nbinning;
00908       gTools().ReadAttr( binch, "NBins", nbinning );
00909       TVectorD binns(nbinning+1);
00910       if (nbinning != nbins) {
00911          Log() << kFATAL << "Number of bins in content and binning array differs"<<Endl;
00912       }
00913       const char* binString = gTools().GetContent(binch);
00914       std::stringstream sb(binString);
00915       for (UInt_t i=0; i<=nbins; i++) sb >> binns[i];
00916       newhist =  new TH1F( hname, hname, nbins, binns.GetMatrixArray() );
00917       newhist->SetDirectory(0);
00918       for (UInt_t i=0; i<nbins; i++) {
00919          s >> val;
00920          newhist->SetBinContent(i+1,val);
00921       }
00922    }
00923 
00924    TString hnameSmooth = hname;
00925    hnameSmooth.ReplaceAll( "_original", "_smoothed" );
00926 
00927    if (fHistOriginal != 0) delete fHistOriginal;
00928    fHistOriginal = newhist;
00929    fHist = (TH1F*)fHistOriginal->Clone( hnameSmooth );
00930    fHist->SetTitle( hnameSmooth );
00931    fHist->SetDirectory(0);
00932 
00933    if (fInterpolMethod == PDF::kKDE) BuildKDEPDF();
00934    else                              BuildSplinePDF();
00935 }
00936 
00937 //_______________________________________________________________________
00938 ostream& TMVA::operator<< ( ostream& os, const PDF& pdf )
00939 {
00940    // write the pdf
00941    Int_t dp = os.precision();
00942    os << "MinNSmooth      " << pdf.fMinNsmooth << std::endl;
00943    os << "MaxNSmooth      " << pdf.fMaxNsmooth << std::endl;
00944    os << "InterpolMethod  " << pdf.fInterpolMethod << std::endl;
00945    os << "KDE_type        " << pdf.fKDEtype << std::endl;
00946    os << "KDE_iter        " << pdf.fKDEiter << std::endl;
00947    os << "KDE_border      " << pdf.fKDEborder << std::endl;
00948    os << "KDE_finefactor  " << pdf.fFineFactor << std::endl;
00949 
00950    TH1* histToWrite = pdf.GetOriginalHist();
00951 
00952    const Int_t nBins = histToWrite->GetNbinsX();
00953 
00954    // note: this is a schema change introduced for v3.7.3
00955    os << "Histogram       "
00956       << histToWrite->GetName()
00957       << "   " << nBins                                 // nbins
00958       << "   " << std::setprecision(12) << histToWrite->GetXaxis()->GetXmin()    // x_min
00959       << "   " << std::setprecision(12) << histToWrite->GetXaxis()->GetXmax()    // x_max
00960       << std::endl;
00961 
00962    // write the smoothed hist
00963    os << "Weights " << std::endl;
00964    os << std::setprecision(8);
00965    for (Int_t i=0; i<nBins; i++) {
00966       os << std::setw(15) << std::left << histToWrite->GetBinContent(i+1) << std::right << " ";
00967       if ((i+1)%5==0) os << std::endl;
00968    }
00969 
00970    os << std::setprecision(dp);
00971    return os; // Return the output stream.
00972 }
00973 
00974 //_______________________________________________________________________
00975 istream& TMVA::operator>> ( istream& istr, PDF& pdf )
00976 {
00977    // read the tree from an istream
00978    TString devnullS;
00979    Int_t   valI;
00980    Int_t   nbins;
00981    Float_t xmin, xmax;
00982    TString hname="_original";
00983    Bool_t doneReading = kFALSE;
00984    while (!doneReading) {
00985       istr >> devnullS;
00986       if (devnullS=="NSmooth")
00987          {istr >> pdf.fMinNsmooth; pdf.fMaxNsmooth=pdf.fMinNsmooth;}
00988       else if (devnullS=="MinNSmooth") istr >> pdf.fMinNsmooth;
00989       else if (devnullS=="MaxNSmooth") istr >> pdf.fMaxNsmooth;
00990       // have to do this with strings to be more stable with developing code
00991       else if (devnullS == "InterpolMethod") { istr >> valI; pdf.fInterpolMethod = PDF::EInterpolateMethod(valI);}
00992       else if (devnullS == "KDE_type")       { istr >> valI; pdf.fKDEtype        = KDEKernel::EKernelType(valI); }
00993       else if (devnullS == "KDE_iter")       { istr >> valI; pdf.fKDEiter        = KDEKernel::EKernelIter(valI);}
00994       else if (devnullS == "KDE_border")     { istr >> valI; pdf.fKDEborder      = KDEKernel::EKernelBorder(valI);}
00995       else if (devnullS == "KDE_finefactor") {
00996          istr  >> pdf.fFineFactor;
00997          if (pdf.GetReadingVersion() != 0 && pdf.GetReadingVersion() < TMVA_VERSION(3,7,3)) { // here we expect the histogram limits if the version is below 3.7.3. When version == 0, the newest TMVA version is assumed.
00998             // coverity[tainted_data_argument]
00999             istr  >> nbins >> xmin >> xmax;
01000             doneReading = kTRUE;
01001          }
01002       }
01003       else if (devnullS == "Histogram")     { istr  >> hname >> nbins >> xmin >> xmax; }
01004       else if (devnullS == "Weights")       { doneReading = kTRUE; }
01005    }
01006 
01007    TString hnameSmooth = hname;
01008    hnameSmooth.ReplaceAll( "_original", "_smoothed" );
01009 
01010    // recreate the original hist
01011    TH1* newhist = new TH1F( hname,hname, nbins, xmin, xmax );
01012    newhist->SetDirectory(0);
01013    Float_t val;
01014    for (Int_t i=0; i<nbins; i++) {
01015       istr >> val;
01016       newhist->SetBinContent(i+1,val);
01017    }
01018 
01019    if (pdf.fHistOriginal != 0) delete pdf.fHistOriginal;
01020    pdf.fHistOriginal = newhist;
01021    pdf.fHist = (TH1F*)pdf.fHistOriginal->Clone( hnameSmooth );
01022    pdf.fHist->SetTitle( hnameSmooth );
01023    pdf.fHist->SetDirectory(0);
01024 
01025    if (pdf.fMinNsmooth>=0) pdf.BuildSplinePDF();
01026    else {
01027       pdf.fInterpolMethod = PDF::kKDE;
01028       pdf.BuildKDEPDF();
01029    }
01030 
01031    return istr;
01032 }
01033 
01034 TMVA::PDF*  TMVA::PDF::ThisPDF( void )
01035 {
01036    // return global "this" pointer of PDF
01037    return fgThisPDF;
01038 }
Generated on Tue Jul 5 15:25:29 2011 for ROOT_528-00b_version by  doxygen 1.5.1