#include <algorithm>
#include <cstdlib>
#include <errno.h>

#include "TObjString.h"
#include "TMath.h"
#include "TString.h"
#include "TTree.h"
#include "TLeaf.h"
#include "TH1.h"
#include "TH2.h"
#include "TList.h"
#include "TSpline.h"
#include "TVector.h"
#include "TMatrixD.h"
#include "TMatrixDSymEigen.h"
#include "TVectorD.h"
#include "TTreeFormula.h"
#include "TXMLEngine.h"
#include "TROOT.h"
#include "TMatrixDSymEigen.h"
#include "TColor.h"
#include "TMVA/Config.h"


#ifndef ROOT_TMVA_Tools
#include "TMVA/Tools.h"
#endif
#ifndef ROOT_TMVA_ROCCalc
#include "TMVA/ROCCalc.h"
#endif
#ifndef ROOT_TMVA_Config
#include "TMVA/Config.h"
#endif
#ifndef ROOT_TMVA_Event
#include "TMVA/Event.h"
#endif
#ifndef ROOT_TMVA_Version
#include "TMVA/Version.h"
#endif
#ifndef ROOT_TMVA_PDF
#include "TMVA/PDF.h"
#endif
#ifndef ROOT_TMVA_MsgLogger
#include "TMVA/MsgLogger.h"
#endif

#include "TMVA/PDF.h"
#include "TMVA/TSpline1.h"
#include "TMVA/TSpline2.h"

using namespace std;

//_______________________________________________________________________________________
TMVA::ROCCalc::ROCCalc(TH1* mvaS, TH1* mvaB) :
   fMaxIter(100),
   fAbsTol(0.0),
   fmvaS(0),
   fmvaB(0),
   fmvaSpdf(0),
   fmvaBpdf(0),
   fSplS(0),
   fSplB(0),
   fSplmvaCumS(0),
   fSplmvaCumB(0),
   fSpleffBvsS(0),
   fnStot(0),
   fnBtot(0),
   fSignificance(0),
   fPurity(0),
   fLogger ( new TMVA::MsgLogger("ROCCalc") )
{
   fUseSplines = kTRUE;
   fNbins      = 100;
   // fmvaS = (TH1*) mvaS->Clone("MVA Signal"); fmvaS->SetTitle("MVA Signal"); 
   // fmvaB = (TH1*) mvaB->Clone("MVA Backgr"); fmvaB->SetTitle("MVA Backgr");
   fmvaS =  mvaS; fmvaS->SetTitle("MVA Signal");
   fmvaB =  mvaB; fmvaB->SetTitle("MVA Backgr");
   fXmax = fmvaS->GetXaxis()->GetXmax();
   fXmin = fmvaS->GetXaxis()->GetXmin(); 
 
   if (TMath::Abs(fXmax-fmvaB->GetXaxis()->GetXmax()) > 0.000001 || 
       TMath::Abs(fXmin-fmvaB->GetXaxis()->GetXmin()) > 0.000001 || 
       fmvaB->GetNbinsX() != fmvaS->GetNbinsX()) {
      Log() << kFATAL << " Cannot cal ROC curve etc, as in put mvaS and mvaB have differen #nbins or range "<<Endl;
   }
   if (!strcmp(fmvaS->GetXaxis()->GetTitle(),"")) fmvaS->SetXTitle("MVA-value");
   if (!strcmp(fmvaB->GetXaxis()->GetTitle(),"")) fmvaB->SetXTitle("MVA-value");
   if (!strcmp(fmvaS->GetYaxis()->GetTitle(),"")) fmvaS->SetYTitle("#entries");
   if (!strcmp(fmvaB->GetYaxis()->GetTitle(),"")) fmvaB->SetYTitle("#entries");
   ApplySignalAndBackgroundStyle(fmvaS, fmvaB);
   fmvaSpdf = mvaS->RebinX(mvaS->GetNbinsX()/100,"MVA Signal PDF"); 
   fmvaBpdf = mvaB->RebinX(mvaB->GetNbinsX()/100,"MVA Backgr PDF");
   fmvaSpdf->SetTitle("MVA Signal PDF"); 
   fmvaBpdf->SetTitle("MVA Backgr PDF");
   fmvaSpdf->Scale(1./fmvaSpdf->GetSumOfWeights());
   fmvaBpdf->Scale(1./fmvaBpdf->GetSumOfWeights());
   fmvaSpdf->SetMaximum(TMath::Max(fmvaSpdf->GetMaximum(), fmvaBpdf->GetMaximum()));
   fmvaBpdf->SetMaximum(TMath::Max(fmvaSpdf->GetMaximum(), fmvaBpdf->GetMaximum()));
   ApplySignalAndBackgroundStyle(fmvaSpdf, fmvaBpdf);

   fCutOrientation = (fmvaS->GetMean() > fmvaB->GetMean()) ? +1 : -1;

   fNevtS = 0;
  
}

//_________________________________________________________________________
void TMVA::ROCCalc::ApplySignalAndBackgroundStyle( TH1* sig, TH1* bkg, TH1* any ) {
   //  Int_t c_Canvas         = TColor::GetColor( "#f0f0f0" );
   //  Int_t c_FrameFill      = TColor::GetColor( "#fffffd" );
   //  Int_t c_TitleBox       = TColor::GetColor( "#5D6B7D" );
   //  Int_t c_TitleBorder    = TColor::GetColor( "#7D8B9D" );
   //  Int_t c_TitleText      = TColor::GetColor( "#FFFFFF" );
   Int_t c_SignalLine     = TColor::GetColor( "#0000ee" );
   Int_t c_SignalFill     = TColor::GetColor( "#7d99d1" );
   Int_t c_BackgroundLine = TColor::GetColor( "#ff0000" );
   Int_t c_BackgroundFill = TColor::GetColor( "#ff0000" );
   //  Int_t c_NovelBlue      = TColor::GetColor( "#2244a5" );

   //signal
   // const Int_t FillColor__S = 38 + 150; // change of Color Scheme in ROOT-5.16.
   // convince yourself with gROOT->GetListOfColors()->Print()
   Int_t FillColor__S = c_SignalFill;
   Int_t FillStyle__S = 1001;
   Int_t LineColor__S = c_SignalLine;
   Int_t LineWidth__S = 2;

   // background
   //Int_t icolor = gConfig().fVariablePlotting.fUsePaperStyle ? 2 + 100 : 2;
   Int_t FillColor__B = c_BackgroundFill;
   Int_t FillStyle__B = 3554;
   Int_t LineColor__B = c_BackgroundLine;
   Int_t LineWidth__B = 2;

   if (sig != NULL) {
      sig->SetLineColor( LineColor__S );
      sig->SetLineWidth( LineWidth__S );
      sig->SetFillStyle( FillStyle__S );
      sig->SetFillColor( FillColor__S );
   }
 
   if (bkg != NULL) {
      bkg->SetLineColor( LineColor__B );
      bkg->SetLineWidth( LineWidth__B );
      bkg->SetFillStyle( FillStyle__B );
      bkg->SetFillColor( FillColor__B );
   }

   if (any != NULL) {
      any->SetLineColor( LineColor__S );
      any->SetLineWidth( LineWidth__S );
      any->SetFillStyle( FillStyle__S );
      any->SetFillColor( FillColor__S );
   }
}


//_______________________________________________________________________________________
TMVA::ROCCalc::~ROCCalc() {
   // destructor
   
   // delete Splines and all histograms that were created only for internal use
   if (fSplS)            { delete fSplS; fSplS = 0; }
   if (fSplB)            { delete fSplB; fSplB = 0; }
   if (fSpleffBvsS)      { delete fSpleffBvsS; fSpleffBvsS = 0; }
   if (fSplmvaCumS)      { delete fSplmvaCumS; fSplmvaCumS = 0; }
   if (fSplmvaCumB)      { delete fSplmvaCumB; fSplmvaCumB = 0; }
   if (fmvaScumul)       { delete fmvaScumul; }
   if (fmvaBcumul)       { delete fmvaBcumul; }

   delete fLogger;
}

//_______________________________________________________________________________________
TH1D* TMVA::ROCCalc::GetROC(){
   // get the ROC curve

   // first get the cumulative distributions of the mva distribution 
   // --> efficiencies vs cut value
   fNevtS = fmvaS->GetSumOfWeights(); // needed to get the error on the eff.. will only be correct if the histogram is not scaled to "integral == 1" Yet;
   if (fNevtS < 2) {
      Log() << kWARNING << "I guess the mva distributions fed into ROCCalc were already normalized, therefore the calculated error on the efficiency will be incorrect !! " << Endl;
      fNevtS = 0;  // reset to zero --> no error will be calculated on the efficiencies
   }
   fmvaScumul = gTools().GetCumulativeDist(fmvaS);
   fmvaBcumul = gTools().GetCumulativeDist(fmvaB);
   fmvaScumul->Scale( 1.0/TMath::Max(std::numeric_limits<double>::epsilon(),fmvaScumul->GetMaximum()) );
   fmvaBcumul->Scale( 1.0/TMath::Max(std::numeric_limits<double>::epsilon(),fmvaBcumul->GetMaximum()) );
   fmvaScumul->SetMinimum(0);
   fmvaBcumul->SetMinimum(0);
   //   fmvaScumul->Draw("hist");
   //   fmvaBcumul->Draw("histsame");

   // background efficiency versus signal efficiency
   TH1D* effBvsS = new TH1D("effBvsS", "ROC-Curve", fNbins, 0, 1 );
   effBvsS->SetXTitle( "Signal eff" );
   effBvsS->SetYTitle( "Backgr eff" );

   // background rejection (=1-eff.) versus signal efficiency
   TH1D* rejBvsS = new TH1D( "rejBvsS", "ROC-Curve", fNbins, 0, 1 );
   rejBvsS->SetXTitle( "Signal eff" );
   rejBvsS->SetYTitle( "Backgr rejection (1-eff)" );
   
   // inverse background eff (1/eff.) versus signal efficiency
   TH1D* inveffBvsS = new TH1D("invBeffvsSeff", "ROC-Curve" , fNbins, 0, 1 );
   inveffBvsS->SetXTitle( "Signal eff" );
   inveffBvsS->SetYTitle( "Inverse backgr. eff (1/eff)" );

   // use root finder
   // spline background efficiency plot
   // note that there is a bin shift when going from a TH1D object to a TGraph :-(
   if (fUseSplines) {
      fSplmvaCumS  = new TSpline1( "spline2_signal",     new TGraph( fmvaScumul ) );
      fSplmvaCumB  = new TSpline1( "spline2_background", new TGraph( fmvaBcumul ) );
      // verify spline sanity
      gTools().CheckSplines( fmvaScumul, fSplmvaCumS );
      gTools().CheckSplines( fmvaBcumul, fSplmvaCumB );
   }

   Double_t effB = 0;
   for (UInt_t bini=1; bini<=fNbins; bini++) {

      // find cut value corresponding to a given signal efficiency
      Double_t effS = effBvsS->GetBinCenter( bini );
      Double_t cut  = Root( effS );

      // retrieve background efficiency for given cut
      if (fUseSplines) effB = fSplmvaCumB->Eval( cut );
      else             effB = fmvaBcumul->GetBinContent( fmvaBcumul->FindBin( cut ) );

      // and fill histograms
      effBvsS->SetBinContent( bini, effB     );
      rejBvsS->SetBinContent( bini, 1.0-effB );
      if (effB>std::numeric_limits<double>::epsilon())
         inveffBvsS->SetBinContent( bini, 1.0/effB );
   }
   
   // create splines for histogram
   fSpleffBvsS = new TSpline1( "effBvsS", new TGraph( effBvsS ) );
   
   // search for overlap point where, when cutting on it,
   // one would obtain: eff_S = rej_B = 1 - eff_B

   Double_t effS = 0., rejB = 0., effS_ = 0., rejB_ = 0.;
   Int_t    nbins = 5000;
   for (Int_t bini=1; bini<=nbins; bini++) {
     
      // get corresponding signal and background efficiencies
      effS = (bini - 0.5)/Float_t(nbins);
      rejB = 1.0 - fSpleffBvsS->Eval( effS );
     
      // find signal efficiency that corresponds to required background efficiency
      if ((effS - rejB)*(effS_ - rejB_) < 0) break;
      effS_ = effS;
      rejB_ = rejB;
   }
   // find cut that corresponds to signal efficiency and update signal-like criterion
   fSignalCut = Root( 0.5*(effS + effS_) );
   
   return rejBvsS;
}

//_______________________________________________________________________________________
Double_t TMVA::ROCCalc::GetROCIntegral(){
   // code to compute the area under the ROC ( rej-vs-eff ) curve

   Double_t effS = 0, effB = 0;
   Int_t    nbins = 1000;
   if (fSpleffBvsS == 0) this->GetROC(); // that will make the ROC calculation if not done yet

   // compute area of rej-vs-eff plot
   Double_t integral = 0;
   for (Int_t bini=1; bini<=nbins; bini++) {
    
      // get corresponding signal and background efficiencies
      effS = (bini - 0.5)/Float_t(nbins);
      effB = fSpleffBvsS->Eval( effS );
      integral += (1.0 - effB);
   }
   integral /= nbins;
  
   return integral;
}

//_______________________________________________________________________________________
Double_t TMVA::ROCCalc::GetEffSForEffBof(Double_t effBref, Double_t &effSerr){
   // get the signal efficiency for a particular backgroud efficiency 
   // that will be the value of the efficiency retured (does not affect
   // the efficiency-vs-bkg plot which is done anyway.
  
   // find precise efficiency value
   Double_t effS=0., effB, effSOld=1., effBOld=0.;
   Int_t    nbins = 1000;
   if (fSpleffBvsS == 0) this->GetROC(); // that will make the ROC calculation if not done yet

   Float_t step=1./nbins;  // stepsize in efficiency binning
   for (Int_t bini=1; bini<=nbins; bini++) {
      // get corresponding signal and background efficiencies
      effS = (bini - 0.5)*step;  // efficiency goes from 0-to-1 in nbins steps of 1/nbins (take middle of the bin)
      effB = fSpleffBvsS->Eval( effS );

      // find signal efficiency that corresponds to required background efficiency
      if ((effB - effBref)*(effBOld - effBref) <= 0) break;
      effSOld = effS;
      effBOld = effB;
   }
  
   // take mean between bin above and bin below
   effS = 0.5*(effS + effSOld);
  
  
   if (fNevtS > 0) effSerr = TMath::Sqrt( effS*(1.0 - effS)/fNevtS );
   else effSerr = 0;

   return effS;
}

//_______________________________________________________________________________________
Double_t TMVA::ROCCalc::GetEffForRoot( Double_t theCut )
{
   // returns efficiency as function of cut
   Double_t retVal=0;

   // retrieve the class object
   if (fUseSplines) retVal = fSplmvaCumS->Eval( theCut );
   else             retVal = fmvaScumul->GetBinContent( fmvaScumul->FindBin( theCut ) );
   
   // caution: here we take some "forbidden" action to hide a problem:
   // in some cases, in particular for likelihood, the binned efficiency distributions
   // do not equal 1, at xmin, and 0 at xmax; of course, in principle we have the
   // unbinned information available in the trees, but the unbinned minimization is
   // too slow, and we don't need to do a precision measurement here. Hence, we force
   // this property.
   Double_t eps = 1.0e-5;
   if      (theCut-fXmin < eps) retVal = (fCutOrientation > 0) ? 1.0 : 0.0;
   else if (fXmax-theCut < eps) retVal = (fCutOrientation > 0) ? 0.0 : 1.0;


   return retVal;
}

//_______________________________________________________________________
Double_t TMVA::ROCCalc::Root( Double_t refValue  )
{
   // Root finding using Brents algorithm; taken from CERNLIB function RZERO
   Double_t a  = fXmin, b = fXmax;
   Double_t fa = GetEffForRoot( a ) - refValue;
   Double_t fb = GetEffForRoot( b ) - refValue;
   if (fb*fa > 0) {
      Log() << kWARNING << "<ROCCalc::Root> initial interval w/o root: "
            << "(a=" << a << ", b=" << b << "),"
            << " (Eff_a=" << GetEffForRoot( a ) 
            << ", Eff_b=" << GetEffForRoot( b ) << "), "
            << "(fa=" << fa << ", fb=" << fb << "), "
            << "refValue = " << refValue << Endl;
      return 1;
   }

   Bool_t   ac_equal(kFALSE);
   Double_t fc = fb;
   Double_t c  = 0, d = 0, e = 0;
   for (Int_t iter= 0; iter <= fMaxIter; iter++) {
      if ((fb < 0 && fc < 0) || (fb > 0 && fc > 0)) {

         // Rename a,b,c and adjust bounding interval d
         ac_equal = kTRUE;
         c  = a; fc = fa;
         d  = b - a; e  = b - a;
      }
  
      if (TMath::Abs(fc) < TMath::Abs(fb)) {
         ac_equal = kTRUE;
         a  = b;  b  = c;  c  = a;
         fa = fb; fb = fc; fc = fa;
      }

      Double_t tol = 0.5 * 2.2204460492503131e-16 * TMath::Abs(b);
      Double_t m   = 0.5 * (c - b);
      if (fb == 0 || TMath::Abs(m) <= tol || TMath::Abs(fb) < fAbsTol) return b;
  
      // Bounds decreasing too slowly: use bisection
      if (TMath::Abs (e) < tol || TMath::Abs (fa) <= TMath::Abs (fb)) { d = m; e = m; }      
      else {
         // Attempt inverse cubic interpolation
         Double_t p, q, r;
         Double_t s = fb / fa;
      
         if (ac_equal) { p = 2 * m * s; q = 1 - s; }
         else {
            q = fa / fc; r = fb / fc;
            p = s * (2 * m * q * (q - r) - (b - a) * (r - 1));
            q = (q - 1) * (r - 1) * (s - 1);
         }
         // Check whether we are in bounds
         if (p > 0) q = -q;
         else       p = -p;
      
         Double_t min1 = 3 * m * q - TMath::Abs (tol * q);
         Double_t min2 = TMath::Abs (e * q);
         if (2 * p < (min1 < min2 ? min1 : min2)) {
            // Accept the interpolation
            e = d;        d = p / q;
         }
         else { d = m; e = m; } // Interpolation failed: use bisection.
      }
      // Move last best guess to a
      a  = b; fa = fb;
      // Evaluate new trial root
      if (TMath::Abs(d) > tol) b += d;
      else                     b += (m > 0 ? +tol : -tol);

      fb = GetEffForRoot( b ) - refValue;

   }

   // Return our best guess if we run out of iterations
   Log() << kWARNING << "<ROCCalc::Root> maximum iterations (" << fMaxIter 
         << ") reached before convergence" << Endl;

   return b;
}

//_______________________________________________________________________
TH1* TMVA::ROCCalc::GetPurity( Int_t nStot, Int_t nBtot)
{
   if (fnStot!=nStot || fnBtot!=nBtot || !fSignificance) {
      GetSignificance(nStot, nBtot); 
      fnStot=nStot; 
      fnBtot=nBtot; 
   }
   return fPurity;
}
//_______________________________________________________________________
TH1* TMVA::ROCCalc::GetSignificance( Int_t nStot, Int_t nBtot)
{
   if (fnStot==nStot && fnBtot==nBtot && !fSignificance) return fSignificance;
   fnStot=nStot; fnBtot=nBtot;

   fSignificance = (TH1*) fmvaScumul->Clone("Significance"); fSignificance->SetTitle("Significance");
   fSignificance->Reset(); fSignificance->SetFillStyle(0);
   fSignificance->SetXTitle("mva cut value");
   fSignificance->SetYTitle("Stat. significance S/Sqrt(S+B)");
   fSignificance->SetLineColor(2);
   fSignificance->SetLineWidth(5);

   fPurity = (TH1*) fmvaScumul->Clone("Purity"); fPurity->SetTitle("Purity");
   fPurity->Reset(); fPurity->SetFillStyle(0);
   fPurity->SetXTitle("mva cut value");
   fPurity->SetYTitle("Purity: S/(S+B)");
   fPurity->SetLineColor(3);
   fPurity->SetLineWidth(5);
   
   Double_t maxSig=0;
   for (Int_t i=1; i<=fSignificance->GetNbinsX(); i++) {
      Double_t S = fmvaScumul->GetBinContent( i ) * nStot;
      Double_t B = fmvaBcumul->GetBinContent( i ) * nBtot;
      Double_t purity;
      Double_t sig;
      if (S+B > 0){
         purity = S/(S+B);
         sig = S/TMath::Sqrt(S+B);
         if (sig > maxSig) {
            maxSig    = sig;
         }
      } else {
         purity=0;
         sig=0;
      }
      cout << "S="<<S<<" B="<<B<< " purity="<<purity<< endl;
      fPurity->SetBinContent( i, purity );
      fSignificance->SetBinContent( i, sig );
   }

   /*   
        TLatex*  line1;
        TLatex*  line2;
        TLatex tl;
        tl.SetNDC();
        tl.SetTextSize( 0.033 );
        Int_t maxbin = fSignificance->GetMaximumBin();
        line1 = tl.DrawLatex( 0.15, 0.23, Form("For %1.0f signal and %1.0f background", nStot, nBtot));
        tl.DrawLatex( 0.15, 0.19, "events the maximum S/Sqrt(S+B) is");

        line2 = tl.DrawLatex( 0.15, 0.15, Form("%4.2f when cutting at %5.2f", 
        maxSig, 
        fSignificance->GetXaxis()->GetBinCenter(maxbin)) );
   */   
   return fSignificance;

}








 ROCCalc.cxx:1
 ROCCalc.cxx:2
 ROCCalc.cxx:3
 ROCCalc.cxx:4
 ROCCalc.cxx:5
 ROCCalc.cxx:6
 ROCCalc.cxx:7
 ROCCalc.cxx:8
 ROCCalc.cxx:9
 ROCCalc.cxx:10
 ROCCalc.cxx:11
 ROCCalc.cxx:12
 ROCCalc.cxx:13
 ROCCalc.cxx:14
 ROCCalc.cxx:15
 ROCCalc.cxx:16
 ROCCalc.cxx:17
 ROCCalc.cxx:18
 ROCCalc.cxx:19
 ROCCalc.cxx:20
 ROCCalc.cxx:21
 ROCCalc.cxx:22
 ROCCalc.cxx:23
 ROCCalc.cxx:24
 ROCCalc.cxx:25
 ROCCalc.cxx:26
 ROCCalc.cxx:27
 ROCCalc.cxx:28
 ROCCalc.cxx:29
 ROCCalc.cxx:30
 ROCCalc.cxx:31
 ROCCalc.cxx:32
 ROCCalc.cxx:33
 ROCCalc.cxx:34
 ROCCalc.cxx:35
 ROCCalc.cxx:36
 ROCCalc.cxx:37
 ROCCalc.cxx:38
 ROCCalc.cxx:39
 ROCCalc.cxx:40
 ROCCalc.cxx:41
 ROCCalc.cxx:42
 ROCCalc.cxx:43
 ROCCalc.cxx:44
 ROCCalc.cxx:45
 ROCCalc.cxx:46
 ROCCalc.cxx:47
 ROCCalc.cxx:48
 ROCCalc.cxx:49
 ROCCalc.cxx:50
 ROCCalc.cxx:51
 ROCCalc.cxx:52
 ROCCalc.cxx:53
 ROCCalc.cxx:54
 ROCCalc.cxx:55
 ROCCalc.cxx:56
 ROCCalc.cxx:57
 ROCCalc.cxx:58
 ROCCalc.cxx:59
 ROCCalc.cxx:60
 ROCCalc.cxx:61
 ROCCalc.cxx:62
 ROCCalc.cxx:63
 ROCCalc.cxx:64
 ROCCalc.cxx:65
 ROCCalc.cxx:66
 ROCCalc.cxx:67
 ROCCalc.cxx:68
 ROCCalc.cxx:69
 ROCCalc.cxx:70
 ROCCalc.cxx:71
 ROCCalc.cxx:72
 ROCCalc.cxx:73
 ROCCalc.cxx:74
 ROCCalc.cxx:75
 ROCCalc.cxx:76
 ROCCalc.cxx:77
 ROCCalc.cxx:78
 ROCCalc.cxx:79
 ROCCalc.cxx:80
 ROCCalc.cxx:81
 ROCCalc.cxx:82
 ROCCalc.cxx:83
 ROCCalc.cxx:84
 ROCCalc.cxx:85
 ROCCalc.cxx:86
 ROCCalc.cxx:87
 ROCCalc.cxx:88
 ROCCalc.cxx:89
 ROCCalc.cxx:90
 ROCCalc.cxx:91
 ROCCalc.cxx:92
 ROCCalc.cxx:93
 ROCCalc.cxx:94
 ROCCalc.cxx:95
 ROCCalc.cxx:96
 ROCCalc.cxx:97
 ROCCalc.cxx:98
 ROCCalc.cxx:99
 ROCCalc.cxx:100
 ROCCalc.cxx:101
 ROCCalc.cxx:102
 ROCCalc.cxx:103
 ROCCalc.cxx:104
 ROCCalc.cxx:105
 ROCCalc.cxx:106
 ROCCalc.cxx:107
 ROCCalc.cxx:108
 ROCCalc.cxx:109
 ROCCalc.cxx:110
 ROCCalc.cxx:111
 ROCCalc.cxx:112
 ROCCalc.cxx:113
 ROCCalc.cxx:114
 ROCCalc.cxx:115
 ROCCalc.cxx:116
 ROCCalc.cxx:117
 ROCCalc.cxx:118
 ROCCalc.cxx:119
 ROCCalc.cxx:120
 ROCCalc.cxx:121
 ROCCalc.cxx:122
 ROCCalc.cxx:123
 ROCCalc.cxx:124
 ROCCalc.cxx:125
 ROCCalc.cxx:126
 ROCCalc.cxx:127
 ROCCalc.cxx:128
 ROCCalc.cxx:129
 ROCCalc.cxx:130
 ROCCalc.cxx:131
 ROCCalc.cxx:132
 ROCCalc.cxx:133
 ROCCalc.cxx:134
 ROCCalc.cxx:135
 ROCCalc.cxx:136
 ROCCalc.cxx:137
 ROCCalc.cxx:138
 ROCCalc.cxx:139
 ROCCalc.cxx:140
 ROCCalc.cxx:141
 ROCCalc.cxx:142
 ROCCalc.cxx:143
 ROCCalc.cxx:144
 ROCCalc.cxx:145
 ROCCalc.cxx:146
 ROCCalc.cxx:147
 ROCCalc.cxx:148
 ROCCalc.cxx:149
 ROCCalc.cxx:150
 ROCCalc.cxx:151
 ROCCalc.cxx:152
 ROCCalc.cxx:153
 ROCCalc.cxx:154
 ROCCalc.cxx:155
 ROCCalc.cxx:156
 ROCCalc.cxx:157
 ROCCalc.cxx:158
 ROCCalc.cxx:159
 ROCCalc.cxx:160
 ROCCalc.cxx:161
 ROCCalc.cxx:162
 ROCCalc.cxx:163
 ROCCalc.cxx:164
 ROCCalc.cxx:165
 ROCCalc.cxx:166
 ROCCalc.cxx:167
 ROCCalc.cxx:168
 ROCCalc.cxx:169
 ROCCalc.cxx:170
 ROCCalc.cxx:171
 ROCCalc.cxx:172
 ROCCalc.cxx:173
 ROCCalc.cxx:174
 ROCCalc.cxx:175
 ROCCalc.cxx:176
 ROCCalc.cxx:177
 ROCCalc.cxx:178
 ROCCalc.cxx:179
 ROCCalc.cxx:180
 ROCCalc.cxx:181
 ROCCalc.cxx:182
 ROCCalc.cxx:183
 ROCCalc.cxx:184
 ROCCalc.cxx:185
 ROCCalc.cxx:186
 ROCCalc.cxx:187
 ROCCalc.cxx:188
 ROCCalc.cxx:189
 ROCCalc.cxx:190
 ROCCalc.cxx:191
 ROCCalc.cxx:192
 ROCCalc.cxx:193
 ROCCalc.cxx:194
 ROCCalc.cxx:195
 ROCCalc.cxx:196
 ROCCalc.cxx:197
 ROCCalc.cxx:198
 ROCCalc.cxx:199
 ROCCalc.cxx:200
 ROCCalc.cxx:201
 ROCCalc.cxx:202
 ROCCalc.cxx:203
 ROCCalc.cxx:204
 ROCCalc.cxx:205
 ROCCalc.cxx:206
 ROCCalc.cxx:207
 ROCCalc.cxx:208
 ROCCalc.cxx:209
 ROCCalc.cxx:210
 ROCCalc.cxx:211
 ROCCalc.cxx:212
 ROCCalc.cxx:213
 ROCCalc.cxx:214
 ROCCalc.cxx:215
 ROCCalc.cxx:216
 ROCCalc.cxx:217
 ROCCalc.cxx:218
 ROCCalc.cxx:219
 ROCCalc.cxx:220
 ROCCalc.cxx:221
 ROCCalc.cxx:222
 ROCCalc.cxx:223
 ROCCalc.cxx:224
 ROCCalc.cxx:225
 ROCCalc.cxx:226
 ROCCalc.cxx:227
 ROCCalc.cxx:228
 ROCCalc.cxx:229
 ROCCalc.cxx:230
 ROCCalc.cxx:231
 ROCCalc.cxx:232
 ROCCalc.cxx:233
 ROCCalc.cxx:234
 ROCCalc.cxx:235
 ROCCalc.cxx:236
 ROCCalc.cxx:237
 ROCCalc.cxx:238
 ROCCalc.cxx:239
 ROCCalc.cxx:240
 ROCCalc.cxx:241
 ROCCalc.cxx:242
 ROCCalc.cxx:243
 ROCCalc.cxx:244
 ROCCalc.cxx:245
 ROCCalc.cxx:246
 ROCCalc.cxx:247
 ROCCalc.cxx:248
 ROCCalc.cxx:249
 ROCCalc.cxx:250
 ROCCalc.cxx:251
 ROCCalc.cxx:252
 ROCCalc.cxx:253
 ROCCalc.cxx:254
 ROCCalc.cxx:255
 ROCCalc.cxx:256
 ROCCalc.cxx:257
 ROCCalc.cxx:258
 ROCCalc.cxx:259
 ROCCalc.cxx:260
 ROCCalc.cxx:261
 ROCCalc.cxx:262
 ROCCalc.cxx:263
 ROCCalc.cxx:264
 ROCCalc.cxx:265
 ROCCalc.cxx:266
 ROCCalc.cxx:267
 ROCCalc.cxx:268
 ROCCalc.cxx:269
 ROCCalc.cxx:270
 ROCCalc.cxx:271
 ROCCalc.cxx:272
 ROCCalc.cxx:273
 ROCCalc.cxx:274
 ROCCalc.cxx:275
 ROCCalc.cxx:276
 ROCCalc.cxx:277
 ROCCalc.cxx:278
 ROCCalc.cxx:279
 ROCCalc.cxx:280
 ROCCalc.cxx:281
 ROCCalc.cxx:282
 ROCCalc.cxx:283
 ROCCalc.cxx:284
 ROCCalc.cxx:285
 ROCCalc.cxx:286
 ROCCalc.cxx:287
 ROCCalc.cxx:288
 ROCCalc.cxx:289
 ROCCalc.cxx:290
 ROCCalc.cxx:291
 ROCCalc.cxx:292
 ROCCalc.cxx:293
 ROCCalc.cxx:294
 ROCCalc.cxx:295
 ROCCalc.cxx:296
 ROCCalc.cxx:297
 ROCCalc.cxx:298
 ROCCalc.cxx:299
 ROCCalc.cxx:300
 ROCCalc.cxx:301
 ROCCalc.cxx:302
 ROCCalc.cxx:303
 ROCCalc.cxx:304
 ROCCalc.cxx:305
 ROCCalc.cxx:306
 ROCCalc.cxx:307
 ROCCalc.cxx:308
 ROCCalc.cxx:309
 ROCCalc.cxx:310
 ROCCalc.cxx:311
 ROCCalc.cxx:312
 ROCCalc.cxx:313
 ROCCalc.cxx:314
 ROCCalc.cxx:315
 ROCCalc.cxx:316
 ROCCalc.cxx:317
 ROCCalc.cxx:318
 ROCCalc.cxx:319
 ROCCalc.cxx:320
 ROCCalc.cxx:321
 ROCCalc.cxx:322
 ROCCalc.cxx:323
 ROCCalc.cxx:324
 ROCCalc.cxx:325
 ROCCalc.cxx:326
 ROCCalc.cxx:327
 ROCCalc.cxx:328
 ROCCalc.cxx:329
 ROCCalc.cxx:330
 ROCCalc.cxx:331
 ROCCalc.cxx:332
 ROCCalc.cxx:333
 ROCCalc.cxx:334
 ROCCalc.cxx:335
 ROCCalc.cxx:336
 ROCCalc.cxx:337
 ROCCalc.cxx:338
 ROCCalc.cxx:339
 ROCCalc.cxx:340
 ROCCalc.cxx:341
 ROCCalc.cxx:342
 ROCCalc.cxx:343
 ROCCalc.cxx:344
 ROCCalc.cxx:345
 ROCCalc.cxx:346
 ROCCalc.cxx:347
 ROCCalc.cxx:348
 ROCCalc.cxx:349
 ROCCalc.cxx:350
 ROCCalc.cxx:351
 ROCCalc.cxx:352
 ROCCalc.cxx:353
 ROCCalc.cxx:354
 ROCCalc.cxx:355
 ROCCalc.cxx:356
 ROCCalc.cxx:357
 ROCCalc.cxx:358
 ROCCalc.cxx:359
 ROCCalc.cxx:360
 ROCCalc.cxx:361
 ROCCalc.cxx:362
 ROCCalc.cxx:363
 ROCCalc.cxx:364
 ROCCalc.cxx:365
 ROCCalc.cxx:366
 ROCCalc.cxx:367
 ROCCalc.cxx:368
 ROCCalc.cxx:369
 ROCCalc.cxx:370
 ROCCalc.cxx:371
 ROCCalc.cxx:372
 ROCCalc.cxx:373
 ROCCalc.cxx:374
 ROCCalc.cxx:375
 ROCCalc.cxx:376
 ROCCalc.cxx:377
 ROCCalc.cxx:378
 ROCCalc.cxx:379
 ROCCalc.cxx:380
 ROCCalc.cxx:381
 ROCCalc.cxx:382
 ROCCalc.cxx:383
 ROCCalc.cxx:384
 ROCCalc.cxx:385
 ROCCalc.cxx:386
 ROCCalc.cxx:387
 ROCCalc.cxx:388
 ROCCalc.cxx:389
 ROCCalc.cxx:390
 ROCCalc.cxx:391
 ROCCalc.cxx:392
 ROCCalc.cxx:393
 ROCCalc.cxx:394
 ROCCalc.cxx:395
 ROCCalc.cxx:396
 ROCCalc.cxx:397
 ROCCalc.cxx:398
 ROCCalc.cxx:399
 ROCCalc.cxx:400
 ROCCalc.cxx:401
 ROCCalc.cxx:402
 ROCCalc.cxx:403
 ROCCalc.cxx:404
 ROCCalc.cxx:405
 ROCCalc.cxx:406
 ROCCalc.cxx:407
 ROCCalc.cxx:408
 ROCCalc.cxx:409
 ROCCalc.cxx:410
 ROCCalc.cxx:411
 ROCCalc.cxx:412
 ROCCalc.cxx:413
 ROCCalc.cxx:414
 ROCCalc.cxx:415
 ROCCalc.cxx:416
 ROCCalc.cxx:417
 ROCCalc.cxx:418
 ROCCalc.cxx:419
 ROCCalc.cxx:420
 ROCCalc.cxx:421
 ROCCalc.cxx:422
 ROCCalc.cxx:423
 ROCCalc.cxx:424
 ROCCalc.cxx:425
 ROCCalc.cxx:426
 ROCCalc.cxx:427
 ROCCalc.cxx:428
 ROCCalc.cxx:429
 ROCCalc.cxx:430
 ROCCalc.cxx:431
 ROCCalc.cxx:432
 ROCCalc.cxx:433
 ROCCalc.cxx:434
 ROCCalc.cxx:435
 ROCCalc.cxx:436
 ROCCalc.cxx:437
 ROCCalc.cxx:438
 ROCCalc.cxx:439
 ROCCalc.cxx:440
 ROCCalc.cxx:441
 ROCCalc.cxx:442
 ROCCalc.cxx:443
 ROCCalc.cxx:444
 ROCCalc.cxx:445
 ROCCalc.cxx:446
 ROCCalc.cxx:447
 ROCCalc.cxx:448
 ROCCalc.cxx:449
 ROCCalc.cxx:450
 ROCCalc.cxx:451
 ROCCalc.cxx:452
 ROCCalc.cxx:453
 ROCCalc.cxx:454
 ROCCalc.cxx:455
 ROCCalc.cxx:456
 ROCCalc.cxx:457
 ROCCalc.cxx:458
 ROCCalc.cxx:459
 ROCCalc.cxx:460
 ROCCalc.cxx:461
 ROCCalc.cxx:462
 ROCCalc.cxx:463
 ROCCalc.cxx:464
 ROCCalc.cxx:465
 ROCCalc.cxx:466
 ROCCalc.cxx:467
 ROCCalc.cxx:468
 ROCCalc.cxx:469
 ROCCalc.cxx:470
 ROCCalc.cxx:471
 ROCCalc.cxx:472
 ROCCalc.cxx:473
 ROCCalc.cxx:474
 ROCCalc.cxx:475
 ROCCalc.cxx:476
 ROCCalc.cxx:477
 ROCCalc.cxx:478
 ROCCalc.cxx:479
 ROCCalc.cxx:480
 ROCCalc.cxx:481
 ROCCalc.cxx:482
 ROCCalc.cxx:483
 ROCCalc.cxx:484
 ROCCalc.cxx:485
 ROCCalc.cxx:486
 ROCCalc.cxx:487
 ROCCalc.cxx:488
 ROCCalc.cxx:489
 ROCCalc.cxx:490
 ROCCalc.cxx:491
 ROCCalc.cxx:492
 ROCCalc.cxx:493
 ROCCalc.cxx:494
 ROCCalc.cxx:495
 ROCCalc.cxx:496
 ROCCalc.cxx:497
 ROCCalc.cxx:498
 ROCCalc.cxx:499