// // @(#)root/hist:$Id$
// Authors: Bartolomeu Rabacal    07/2010
/**********************************************************************
 *                                                                    *
 * Copyright (c) 2006 , ROOT MathLib Team                             *
 *                                                                    *
 * For the licensing terms see $ROOTSYS/LICENSE.                      *
 * For the list of contributors see $ROOTSYS/README/CREDITS.          *
 *                                                                    *
 **********************************************************************/
//////////////////////////////////////////////////////////////////////////////
/*
 Kernel Density Estimation class.
 The three main references are (1) "Scott DW, Multivariate Density Estimation. Theory, Practice and Visualization. New York: Wiley",
 (2) "Jann Ben - ETH Zurich, Switzerland -, Univariate kernel density estimation document for KDENS:
 Stata module for univariate kernel density estimation."
 (3) "Hardle W, Muller M, Sperlich S, Werwatz A, Nonparametric and Semiparametric Models. Springer."
 The algorithm is briefly described in (4) "Cranmer KS, Kernel Estimation in High-Energy
 Physics. Computer Physics Communications 136:198-207,2001" - e-Print Archive: hep ex/0011057.
 A binned version is also implemented to address the performance issue due to its data size dependance.
 */


#include <functional>
#include <algorithm>
#include <numeric>
#include <limits>
#include <cassert>

#include "Math/Error.h"
#include "TMath.h"
#include "Math/Functor.h"
#include "Math/Integrator.h"
#include "Math/QuantFuncMathCore.h"
#include "Math/RichardsonDerivator.h"
#include "TGraphErrors.h"
#include "TF1.h"
#include "TH1.h"
#include "TCanvas.h"
#include "TKDE.h"


ClassImp(TKDE)

class TKDE::TKernel {
   TKDE* fKDE;
   UInt_t fNWeights; // Number of kernel weights (bandwidth as vectorized for binning)
   std::vector<Double_t> fWeights; // Kernel weights (bandwidth)
public:
   TKernel(Double_t weight, TKDE* kde);
   void ComputeAdaptiveWeights();
   Double_t operator()(Double_t x) const;
   Double_t GetWeight(Double_t x) const;
   Double_t GetFixedWeight() const;
   const std::vector<Double_t> & GetAdaptiveWeights() const;
};

struct TKDE::KernelIntegrand {
   enum EIntegralResult{kNorm, kMu, kSigma2, kUnitIntegration};
   KernelIntegrand(const TKDE* kde, EIntegralResult intRes);
   Double_t operator()(Double_t x) const;
private:
   const TKDE* fKDE;
   EIntegralResult fIntegralResult;
};

// TKDE::TKDE(UInt_t events, const Double_t* data, Double_t xMin, Double_t xMax, const Option_t* option, Double_t rho) :
//    fData(events, 0.0),
//    fEvents(events, 0.0),
//    fPDF(0),
//    fUpperPDF(0),
//    fLowerPDF(0),
//    fApproximateBias(0),
//    fGraph(0),
//    fNewData(false),
//    fUseMinMaxFromData((xMin >= xMax)),
//    fNBins(events < 10000 ? 100: events / 10),
//    fNEvents(events),
//    fUseBinsNEvents(10000),
//    fMean(0.0),
//    fSigma(0.0),
//    fXMin(xMin),
//    fXMax(xMax),
//    fAdaptiveBandwidthFactor(1.0),
//    fCanonicalBandwidths(std::vector<Double_t>(kTotalKernels, 0.0)),
//    fKernelSigmas2(std::vector<Double_t>(kTotalKernels, -1.0)),
//    fSettedOptions(std::vector<Bool_t>(4, kFALSE))
// {
//    //Class constructor
//    SetOptions(option, rho);
//    CheckOptions();
//    SetMirror();
//    SetUseBins();
//    SetKernelFunction();
//    SetData(data);
//    SetCanonicalBandwidths();
//    SetKernelSigmas2();
//    SetKernel();
// }

TKDE::~TKDE() {
   //Class destructor
   if (fPDF)              delete fPDF;
   if (fUpperPDF)         delete fUpperPDF;
   if (fLowerPDF)         delete fLowerPDF;
   if (fGraph)         delete fGraph;
   if (fApproximateBias)  delete fApproximateBias;
   if (fKernelFunction && fKernelType != kUserDefined) delete fKernelFunction;
   if (fKernel) delete fKernel;
}

void TKDE::Instantiate(KernelFunction_Ptr kernfunc, UInt_t events, const Double_t* data, const Double_t* dataWeights, Double_t xMin, Double_t xMax, const Option_t* option, Double_t rho) {
   // Template's constructor surrogate
   fData = std::vector<Double_t>(events, 0.0);
   fEvents = std::vector<Double_t>(events, 0.0);
   fPDF = 0;
   fKernel = 0;
   fKernelFunction = 0;
   fUpperPDF = 0;
   fLowerPDF = 0;
   fApproximateBias = 0;
   fGraph = 0;
   fNewData = false;
   fUseMirroring = false; fMirrorLeft = false; fMirrorRight = false;
   fAsymLeft = false; fAsymRight = false; 
   fNBins = events < 10000 ? 100 : events / 10;
   fNEvents = events;
   fUseBinsNEvents = 10000;
   fMean = 0.0;
   fSigma = 0.0;
   fXMin = xMin;
   fXMax = xMax;
   fUseMinMaxFromData = (fXMin >= fXMax);
   fSumOfCounts = 0;
   fAdaptiveBandwidthFactor = 1.;
   fRho = rho;
   fWeightSize = 0;
   fCanonicalBandwidths = std::vector<Double_t>(kTotalKernels, 0.0);
   fKernelSigmas2 = std::vector<Double_t>(kTotalKernels, -1.0);
   fSettedOptions = std::vector<Bool_t>(4, kFALSE);
   SetOptions(option, rho);
   CheckOptions(kTRUE);
   SetMirror();
   SetUseBins();
   SetData(data, dataWeights);
   SetKernelFunction(kernfunc);
}

void TKDE::SetOptions(const Option_t* option, Double_t rho) {
   //Sets User defined construction options
   TString opt = option;
   opt.ToLower();
   std::string options = opt.Data();
   size_t numOpt = 4;
   std::vector<std::string> voption(numOpt, "");
   for (std::vector<std::string>::iterator it = voption.begin(); it != voption.end() && !options.empty(); ++it) {
      size_t pos = options.find_last_of(';');
      if (pos == std::string::npos) {
         *it = options;
         break;
      }
      *it = options.substr(pos + 1);
      options = options.substr(0, pos);
   }
   for (std::vector<std::string>::iterator it = voption.begin(); it != voption.end(); ++it) {
      size_t pos = (*it).find(':');
      if (pos != std::string::npos) {
         GetOptions((*it).substr(0, pos) , (*it).substr(pos + 1));
      }
   }
   AssureOptions();
   fRho = rho;
}

void TKDE::SetDrawOptions(const Option_t* option, TString& plotOpt, TString& drawOpt) {
   // Sets User defined drawing options
   size_t numOpt = 2;
   std::string options = TString(option).Data();
   std::vector<std::string> voption(numOpt, "");
   for (std::vector<std::string>::iterator it = voption.begin(); it != voption.end() && !options.empty(); ++it) {
      size_t pos = options.find_last_of(';');
      if (pos == std::string::npos) {
         *it = options;
         break;
      }
      *it = options.substr(pos + 1);
      options = options.substr(0, pos);
   }
   Bool_t foundPlotOPt = kFALSE;
   Bool_t foundDrawOPt = kFALSE;
   for (std::vector<std::string>::iterator it = voption.begin(); it != voption.end() && !options.empty(); ++it) {
      size_t pos = (*it).find(':');
      if (pos == std::string::npos) break;
      TString optionType = (*it).substr(0, pos);
      TString optionInstance = (*it).substr(pos + 1);
      optionType.ToLower();
      optionInstance.ToLower();
      if (optionType.Contains("plot")) {
         foundPlotOPt = kTRUE;
         if (optionInstance.Contains("estimate") || optionInstance.Contains("errors") || optionInstance.Contains("confidenceinterval"))
            plotOpt = optionInstance;
         else
            this->Warning("SetDrawOptions", "Unknown plotting option: setting to KDE estimate plot.");
      } else if (optionType.Contains("drawoptions")) {
         foundDrawOPt = kTRUE;
         drawOpt = optionInstance;
      }
   }
   if (!foundPlotOPt) {
      this->Warning("SetDrawOptions", "No plotting option: setting to KDE estimate plot.");
      plotOpt = "estimate";
   }
   if (!foundDrawOPt) {
      this->Warning("SetDrawOptions", "No drawing options: setting to default ones.");
      drawOpt = "apl4";
   }
}

void TKDE::GetOptions(std::string optionType, std::string option) {
   // Gets User defined KDE construction options
   if (optionType.compare("kerneltype") == 0) {
      fSettedOptions[0] = kTRUE;
      if (option.compare("gaussian") == 0) {
         fKernelType = kGaussian;
      } else if (option.compare("epanechnikov") == 0) {
         fKernelType = kEpanechnikov;
      } else if (option.compare("biweight") == 0) {
         fKernelType = kBiweight;
      } else if (option.compare("cosinearch") == 0) {
         fKernelType = kCosineArch;
      } else if (option.compare("userdefined") == 0) {
         fKernelType = kUserDefined;
      } else {
         this->Warning("GetOptions", "Unknown kernel type option: setting to Gaussian");
         fKernelType = kGaussian;
      }
   } else if (optionType.compare("iteration") == 0) {
      fSettedOptions[1] = kTRUE;
      if (option.compare("adaptive") == 0) {
         fIteration = kAdaptive;
      } else if (option.compare("fixed") == 0) {
         fIteration = kFixed;
      } else {
         this->Warning("GetOptions", "Unknown iteration option: setting to Adaptive");
         fIteration = kAdaptive;
      }
   } else if (optionType.compare("mirror") == 0) {
      fSettedOptions[2] = kTRUE;
      if (option.compare("nomirror") == 0) {
         fMirror = kNoMirror;
      } else if (option.compare("mirrorleft") == 0) {
         fMirror = kMirrorLeft;
      } else if (option.compare("mirrorright") == 0) {
         fMirror = kMirrorRight;
      } else if (option.compare("mirrorboth") == 0) {
         fMirror = kMirrorBoth;
      } else if (option.compare("mirrorasymleft") == 0) {
         fMirror = kMirrorAsymLeft;
      } else if (option.compare("mirrorasymleftright") == 0) {
         fMirror = kMirrorAsymLeftRight;
      } else if (option.compare("mirrorasymright") == 0) {
         fMirror = kMirrorAsymRight;
      } else if (option.compare("mirrorleftasymright") == 0) {
         fMirror = kMirrorLeftAsymRight;
      } else if (option.compare("mirrorasymboth") == 0) {
         fMirror = kMirrorAsymBoth;
      } else {
         this->Warning("GetOptions", "Unknown mirror option: setting to NoMirror");
         fMirror = kNoMirror;
      }
   } else if (optionType.compare("binning") == 0) {
      fSettedOptions[3] = kTRUE;
      if (option.compare("unbinned") == 0) {
         fBinning = kUnbinned;
      } else if (option.compare("relaxedbinning") == 0) {
         fBinning = kRelaxedBinning;
      } else if (option.compare("forcedbinning") == 0) {
         fBinning = kForcedBinning;
      } else {
         this->Warning("GetOptions", "Unknown binning option: setting to RelaxedBinning");
         fBinning = kRelaxedBinning;
      }
   }
}

void TKDE::AssureOptions() {
   // Sets missing construction options to default ones
   if (!fSettedOptions[0]) {
      fKernelType = kGaussian;
   }
   if (!fSettedOptions[1]) {
      fIteration = kAdaptive;
   }
   if (!fSettedOptions[2]) {
      fMirror = kNoMirror;
   }
   if (!fSettedOptions[3]) {
      fBinning = kRelaxedBinning;
   }
}

void TKDE::CheckOptions(Bool_t isUserDefinedKernel) {
   // Sets User global options
   if (!(isUserDefinedKernel) && !(fKernelType >= kGaussian && fKernelType < kUserDefined)) {
      Error("CheckOptions", "Illegal user kernel type input! Use template constructor for user defined kernel.");
   }
   if (fIteration != kAdaptive && fIteration != kFixed) {
      Warning("CheckOptions", "Illegal user iteration type input - use default value !");
      fIteration = kAdaptive;
   }
   if (!(fMirror >= kNoMirror && fMirror <= kMirrorAsymBoth)) {
      Warning("CheckOptions", "Illegal user mirroring type input - use default value !");
      fMirror = kNoMirror;
   }
   if (!(fBinning >= kUnbinned && fBinning <= kForcedBinning)) {
      Warning("CheckOptions", "Illegal user binning type input - use default value !");
      fBinning = kRelaxedBinning;
   }
   if (fRho <= 0.0) {
      Warning("CheckOptions", "Tuning factor rho cannot be non-positive - use default value !");
      fRho = 1.0;
   }
}

void TKDE::SetKernelType(EKernelType kern) {
   // Sets User option for the choice of kernel estimator
   if (fKernelFunction && fKernelType != kUserDefined) {
      delete fKernelFunction;
      fKernelFunction = 0;
   }
   fKernelType = kern;
   CheckOptions();
   SetKernelFunction(0);
}

void TKDE::SetIteration(EIteration iter) {
   // Sets User option for fixed or adaptive iteration
   fIteration = iter;
   CheckOptions();
   SetKernel();
}


void TKDE::SetMirror(EMirror mir) {
   // Sets User option for mirroring the data
   fMirror = mir;
   CheckOptions();
   SetMirror();
   if (fUseMirroring) {
      SetMirroredEvents();
   }
   SetKernel();
}


void TKDE::SetBinning(EBinning bin) {
   // Sets User option for binning the weights
   fBinning = bin;
   CheckOptions();
   SetUseBins();
   SetKernel();
}

void TKDE::SetNBins(UInt_t nbins) {
   // Sets User option for number of bins
   if (!nbins) {
      Error("SetNBins", "Number of bins must be greater than zero.");
      return;
   }
   fNBins = nbins;
   fWeightSize = fNBins / (fXMax - fXMin);
   SetBinCentreData(fXMin, fXMax);
   SetBinCountData();

   if (fBinning == kUnbinned) {
      Warning("SetNBins", "Bin type using SetBinning must set for using a binned evaluation");
      return;
   }

   SetKernel();
}

void TKDE::SetUseBinsNEvents(UInt_t nEvents) {
   // Sets User option for the minimum number of events for allowing automatic binning
   fUseBinsNEvents = nEvents;
   SetUseBins();
   SetKernel();
}

void TKDE::SetTuneFactor(Double_t rho) {
   // Factor which can be used to tune the smoothing.
   // It is used as multiplicative factor for the fixed and adaptive bandwidth.
   // A value < 1 will reproduce better the tails but oversmooth the peak
   // while a factor > 1 will overestimate the tail
   fRho = rho;
   CheckOptions();
   SetKernel();
}

void TKDE::SetRange(Double_t xMin, Double_t xMax) {
   // Sets minimum range value and maximum range value
   if (xMin >= xMax) {
      Error("SetRange", "Minimum range cannot be bigger or equal than the maximum range! Present range values remain the same.");
      return;
   }
   fXMin = xMin;
   fXMax = xMax;
   fUseMinMaxFromData = false;
   SetKernel();
}

// private methods

void TKDE::SetUseBins() {
   // Sets User option for using binned weights
   switch (fBinning) {
      default:
      case kRelaxedBinning:
         if (fNEvents >= fUseBinsNEvents) {
            fUseBins = kTRUE;
         } else {
            fUseBins = kFALSE;
         }
         break;
      case kForcedBinning:
         fUseBins = kTRUE;
         break;
      case kUnbinned:
         fUseBins = kFALSE;
   }
}

void TKDE::SetMirror() {
   // Sets the mirroring
   fMirrorLeft   = fMirror == kMirrorLeft      || fMirror == kMirrorBoth          || fMirror == kMirrorLeftAsymRight;
   fMirrorRight  = fMirror == kMirrorRight     || fMirror == kMirrorBoth          || fMirror == kMirrorAsymLeftRight;
   fAsymLeft     = fMirror == kMirrorAsymLeft  || fMirror == kMirrorAsymLeftRight || fMirror == kMirrorAsymBoth;
   fAsymRight    = fMirror == kMirrorAsymRight || fMirror == kMirrorLeftAsymRight || fMirror == kMirrorAsymBoth;
   fUseMirroring = fMirrorLeft                 || fMirrorRight ;
}

void TKDE::SetData(const Double_t* data, const Double_t* wgts) {
   // Sets the data events input sample or bin centres for binned option and computes basic estimators
   if (!data) {
      if (fNEvents) fData.reserve(fNEvents);
      return;
   }
   fEvents.assign(data, data + fNEvents);
   if (wgts) fEventWeights.assign(wgts, wgts + fNEvents);
   
   if (fUseMinMaxFromData) {
      fXMin = *std::min_element(fEvents.begin(), fEvents.end());
      fXMax = *std::max_element(fEvents.begin(), fEvents.end());
   }

   if (fUseBins) {
      if (fNBins >= fNEvents) {
         this->Warning("SetData", "Default number of bins is greater or equal to number of events. Use SetNBins(UInt_t) to set the appropriate number of bins");
      }
      fWeightSize = fNBins / (fXMax - fXMin);
      SetBinCentreData(fXMin, fXMax);
   } else {
      fWeightSize = fNEvents / (fXMax - fXMin);
      fData = fEvents;
   }
   // to set fBinCOunt and fSumOfCounts
   SetBinCountData();


   ComputeDataStats(); 
   if (fUseMirroring) {
      SetMirroredEvents();
   }
}

void TKDE::InitFromNewData() {
   // re-initialize when new data have been filled in TKDE
   // re-compute kernel quantities and mean and sigma
   fNewData = false;
   fEvents = fData;
   if (fUseMinMaxFromData) {
      fXMin = *std::min_element(fEvents.begin(), fEvents.end());
      fXMax = *std::max_element(fEvents.begin(), fEvents.end());
   }
   ComputeDataStats(); 
   //    if (fUseBins) {
   // } // bin usage is not supported in this case
   //
   fWeightSize = fNEvents / (fXMax - fXMin);
   if (fUseMirroring) {
      SetMirroredEvents();
   }
   SetKernel();
}

void TKDE::SetMirroredEvents() {
   // Mirrors the data 
   std::vector<Double_t> originalEvents = fEvents;
   std::vector<Double_t> originalWeights = fEventWeights;
  if (fMirrorLeft) {
      fEvents.resize(2 * fNEvents, 0.0);
      transform(fEvents.begin(), fEvents.begin() + fNEvents, fEvents.begin() + fNEvents, std::bind1st(std::minus<Double_t>(), 2 * fXMin));
   }
   if (fMirrorRight) {
      fEvents.resize((fMirrorLeft + 2) * fNEvents, 0.0);
      transform(fEvents.begin(), fEvents.begin() + fNEvents, fEvents.begin() + (fMirrorLeft + 1) * fNEvents, std::bind1st(std::minus<Double_t>(), 2 * fXMax));
   }
   if (!fEventWeights.empty() && (fMirrorLeft || fMirrorRight)) {
      // copy weights too
      fEventWeights.insert(fEventWeights.end(), fEventWeights.begin(), fEventWeights.end() ); 
   }

   if(fUseBins) {
      fNBins *= (fMirrorLeft + fMirrorRight + 1);
      Double_t xmin = fMirrorLeft  ? 2 * fXMin - fXMax : fXMin;
      Double_t xmax = fMirrorRight ? 2 * fXMax - fXMin : fXMax;
      SetBinCentreData(xmin, xmax);
   } else {
      fData = fEvents;
   }
   SetBinCountData();
   
   fEvents = originalEvents;
   fEventWeights = originalWeights; 
}

void TKDE::SetMean() {
   // Computes input data's mean
   fMean = std::accumulate(fEvents.begin(), fEvents.end(), 0.0) / fEvents.size();
}

void TKDE::SetSigma(Double_t R) {
   // Computes input data's sigma
   fSigma = std::sqrt(1. / (fEvents.size() - 1.) * (std::inner_product(fEvents.begin(), fEvents.end(), fEvents.begin(), 0.0) - fEvents.size() * std::pow(fMean, 2.)));
   fSigmaRob = std::min(fSigma, R / 1.349); // Sigma's robust estimator
}

void TKDE::SetKernel() {
   // Sets the kernel density estimator
   UInt_t n = fData.size();
   if (n == 0) return;
   // Optimal bandwidth (Silverman's rule of thumb with assumed Gaussian density)
   Double_t weight = fCanonicalBandwidths[kGaussian] * fSigmaRob * std::pow(3. / (8. * std::sqrt(M_PI)) * n, -0.2);
   weight *= fRho * fCanonicalBandwidths[fKernelType] / fCanonicalBandwidths[kGaussian];
   if (fKernel) delete fKernel;
   fKernel = new TKernel(weight, this);
   if (fIteration == kAdaptive) {
      fKernel->ComputeAdaptiveWeights();
   }
}

void TKDE::SetKernelFunction(KernelFunction_Ptr kernfunc) {

   assert(fKernelFunction == 0);  // to avoid memory leaks
   switch (fKernelType) {
      case kGaussian :
         fKernelFunction = new ROOT::Math::WrappedMemFunction<TKDE, Double_t (TKDE::*)(Double_t) const>(*this, &TKDE::GaussianKernel);
         break;
      case kEpanechnikov :
         fKernelFunction = new ROOT::Math::WrappedMemFunction<TKDE, Double_t (TKDE::*)(Double_t) const>(*this, &TKDE::EpanechnikovKernel);
         break;
      case kBiweight :
         fKernelFunction = new ROOT::Math::WrappedMemFunction<TKDE, Double_t (TKDE::*)(Double_t) const>(*this, &TKDE::BiweightKernel);
         break;
      case kCosineArch :
         fKernelFunction = new ROOT::Math::WrappedMemFunction<TKDE, Double_t (TKDE::*)(Double_t) const>(*this, &TKDE::CosineArchKernel);
         break;
      case kUserDefined :
         fKernelFunction = kernfunc;
         if (fKernelFunction)  CheckKernelValidity();
         break;
      case kTotalKernels :
      default:
         /// for user defined kernels
         fKernelFunction = kernfunc;
         fKernelType = kUserDefined; 
   }

   if (fKernelType == kUserDefined) { 
      if (fKernelFunction) {
         CheckKernelValidity();
         SetUserCanonicalBandwidth();
         SetUserKernelSigma2();
      }
      else {
         Error("SetKernelFunction", "User kernel function is not defined !");
         return;
      }
   }
   assert(fKernelFunction); 
   SetKernelSigmas2();
   SetCanonicalBandwidths();
   SetKernel();
}

void TKDE::SetCanonicalBandwidths() {
   // Sets the canonical bandwidths according to the kernel type
   fCanonicalBandwidths[kGaussian] = 0.7764;     // Checked in Mathematica
   fCanonicalBandwidths[kEpanechnikov] = 1.7188; // Checked in Mathematica
   fCanonicalBandwidths[kBiweight] = 2.03617;    // Checked in Mathematica
   fCanonicalBandwidths[kCosineArch] = 1.7663;   // Checked in Mathematica
   fCanonicalBandwidths[kUserDefined] = 1.0;     // To be Checked 
}

void TKDE::SetKernelSigmas2() {
   // Sets the kernel sigmas2 according to the kernel type
   fKernelSigmas2[kGaussian] = 1.0;
   fKernelSigmas2[kEpanechnikov] = 1.0 / 5.0;
   fKernelSigmas2[kBiweight] = 1.0 / 7.0;
   fKernelSigmas2[kCosineArch] = 1.0 - 8.0 / std::pow(M_PI, 2);
}

TF1* TKDE::GetFunction(UInt_t npx, Double_t xMin, Double_t xMax) {
   // Returns the PDF estimate as a function sampled in npx between xMin and xMax
   // the KDE is not re-normalized to the xMin/xMax range.
   // The user manages the returned function
   // For getting a non-sampled TF1, one can create a TF1 directly from the TKDE by doing
   // TF1 * f1  = new TF1("f1",kde,xMin,xMax,0);
   return GetKDEFunction(npx,xMin,xMax);
}

TF1* TKDE::GetUpperFunction(Double_t confidenceLevel, UInt_t npx, Double_t xMin, Double_t xMax) {
   // Returns the PDF upper estimate (upper confidence interval limit)
   return GetPDFUpperConfidenceInterval(confidenceLevel,npx,xMin,xMax);
}

TF1* TKDE::GetLowerFunction(Double_t confidenceLevel, UInt_t npx, Double_t xMin, Double_t xMax) {
   // Returns the PDF lower estimate (lower confidence interval limit)
   return GetPDFLowerConfidenceInterval(confidenceLevel,npx,xMin,xMax);
}

TF1* TKDE::GetApproximateBias(UInt_t npx, Double_t xMin, Double_t xMax) {
   // Returns the PDF estimate bias
   return GetKDEApproximateBias(npx,xMin,xMax);
}

void TKDE::Fill(Double_t data) {
   // Fills data member with User input data event for the unbinned option
   if (fUseBins) {
      this->Warning("Fill", "Cannot fill data with data binned option. Data input ignored.");
      return;
   }
   fData.push_back(data);
   fNEvents++;
   fNewData = kTRUE;
}

void TKDE::Fill(Double_t data, Double_t weight) {
   // Fills data member with User input data event for the unbinned option
   if (fUseBins) {
      this->Warning("Fill", "Cannot fill data with data binned option. Data input ignored.");
      return;
   }
   fData.push_back(data);  // should not be here fEvent ??
   fEventWeights.push_back(weight);
   fNEvents++;
   fNewData = kTRUE;
}

Double_t TKDE::operator()(const Double_t* x, const Double_t*) const {
   // The class's unary function: returns the kernel density estimate
   return (*this)(*x);
}

Double_t TKDE::operator()(Double_t x) const {
   // The class's unary function: returns the kernel density estimate
   if (fNewData) (const_cast<TKDE*>(this))->InitFromNewData();
   return (*fKernel)(x);
}

Double_t TKDE::GetMean() const {
   // return the mean of the data
   if (fNewData) (const_cast<TKDE*>(this))->InitFromNewData();
   return fMean;
}

Double_t TKDE::GetSigma() const {
   // return the standard deviation  of the data
   if (fNewData) (const_cast<TKDE*>(this))->InitFromNewData();
   return fSigma;
}

Double_t TKDE::GetRAMISE() const {
   // Returns the Root Asymptotic Mean Integrated Squared Error according to Silverman's rule of thumb with assumed Gaussian density
   Double_t result = 5. / 4. * fKernelSigmas2[fKernelType] * std::pow(fCanonicalBandwidths[fKernelType], 4) * std::pow(3. / (8. * std::sqrt(M_PI)) , -0.2) * fSigmaRob * std::pow(fNEvents, -0.8);
   return std::sqrt(result);
}

TKDE::TKernel::TKernel(Double_t weight, TKDE* kde) :
// Internal class constructor
fKDE(kde),
fNWeights(kde->fData.size()),
fWeights(fNWeights, weight)
{}

void TKDE::TKernel::ComputeAdaptiveWeights() {
   // Gets the adaptive weights (bandwidths) for TKernel internal computation
   std::vector<Double_t> weights = fWeights;
   Double_t minWeight = weights[0] * 0.05;
   unsigned int n = fKDE->fData.size();
   assert( n == weights.size() );
   bool useDataWeights = (fKDE->fBinCount.size() == n); 
   Double_t f = 0.0;
   for (unsigned int i = 0; i < n; ++i) { 
//   for (; weight != weights.end(); ++weight, ++data, ++dataW) {
      if (useDataWeights && fKDE->fBinCount[i] <= 0) continue;  // skip negative or null weights
      f = (*fKDE->fKernel)(fKDE->fData[i]);
      if (f <= 0)
         fKDE->Warning("ComputeAdativeWeights","function value is zero or negative for x = %f w = %f",
                       fKDE->fData[i],(useDataWeights) ? fKDE->fBinCount[i] : 1.);
      weights[i] = std::max(weights[i] /= std::sqrt(f), minWeight);
      fKDE->fAdaptiveBandwidthFactor += std::log(f);
      //printf("(f = %f w = %f af = %f ),",f,*weight,fKDE->fAdaptiveBandwidthFactor);
   }
   Double_t kAPPROX_GEO_MEAN = 0.241970724519143365; // 1 / TMath::Power(2 * TMath::Pi(), .5) * TMath::Exp(-.5). Approximated geometric mean over pointwise data (the KDE function is substituted by the "real Gaussian" pdf) and proportional to sigma. Used directly when the mirroring is enabled, otherwise computed from the data
   fKDE->fAdaptiveBandwidthFactor = fKDE->fUseMirroring ? kAPPROX_GEO_MEAN / fKDE->fSigmaRob : std::sqrt(std::exp(fKDE->fAdaptiveBandwidthFactor / fKDE->fData.size()));
   transform(weights.begin(), weights.end(), fWeights.begin(), std::bind2nd(std::multiplies<Double_t>(), fKDE->fAdaptiveBandwidthFactor));
   //printf("adaptive bandwidth factor % f weight 0 %f , %f \n",fKDE->fAdaptiveBandwidthFactor, weights[0],fWeights[0] );
}

Double_t TKDE::TKernel::GetWeight(Double_t x) const {
   // Returns the bandwidth
   return fWeights[fKDE->Index(x)];
}

void TKDE::SetBinCentreData(Double_t xmin, Double_t xmax) {
   // Returns the bins' centres from the data for using with the binned option
   fData.assign(fNBins, 0.0);
   Double_t binWidth((xmax - xmin) / fNBins);
   for (UInt_t i = 0; i < fNBins; ++i) {
      fData[i] = xmin + (i + 0.5) * binWidth;
   }
}

void TKDE::SetBinCountData() {
   // Returns the bins' count from the data for using with the binned option
   // or set the bin count to the weights in case of weighted data
   if (fUseBins) { 
      fBinCount.resize(fNBins);
      fSumOfCounts = 0;
      // case of weighted events 
      if (!fEventWeights.empty() ) { 
         for (UInt_t i = 0; i < fNEvents; ++i) {
            if (fEvents[i] >= fXMin && fEvents[i] < fXMax) {
               fBinCount[Index(fEvents[i])] += fEventWeights[i];
               fSumOfCounts += fEventWeights[i];
               //printf("sum of counts %f - bin count %d - %f \n",fSumOfCounts, Index(fEvents[i]), fBinCount[Index(fEvents[i])] );
            }
         }
      }
      // case of unweighted data 
      else {
         for (UInt_t i = 0; i < fNEvents; ++i) {
            if (fEvents[i] >= fXMin && fEvents[i] < fXMax) {
               fBinCount[Index(fEvents[i])] += 1;
               fSumOfCounts += 1;
            }
         }
      }
   }
   else if (!fEventWeights.empty() ) {
      fBinCount = fEventWeights;
      fSumOfCounts = 0;
      for (UInt_t i = 0; i < fNEvents; ++i) 
         fSumOfCounts += fEventWeights[i];
   }
   else {
      fSumOfCounts = fNEvents;
      fBinCount.clear(); 
   }
}

void TKDE::Draw(const Option_t* opt) {
   // Draws either the KDE functions or its errors
   // Possible options:
   //                    ""  (default) - draw just the kde
   //                    "same" draw on top of existing pad
   //                    "Errors" draw a TGraphErrors with the point and errors
   //                    "confidenceinterval" draw KDE + conf interval functions (default is 95%)
   //                    "confidenceinterval@0.90" draw KDE + conf interval functions at 90%
   //                      Extra options can be passed in opt for drawing the TF1 or the TGraph
   //
   //NOTE:  The functions GetDrawnFunction(), GetDrawnUpperFunction(), GetDrawnLowerFunction()
   //  and GetGraphWithErrors() return the corresponding drawn objects (which are maneged by the TKDE)
   // They can be used to changes style, color, etc...

   // TString plotOpt = "";
   // TString drawOpt = "";
   // LM : this is too complicates - skip it - not needed for just
   // three options
   // SetDrawOptions(opt, plotOpt, drawOpt);
   TString plotOpt = opt;
   plotOpt.ToLower();
   TString drawOpt = plotOpt;
   if(gPad && !plotOpt.Contains("same")) {
      gPad->Clear();
   }
   if (plotOpt.Contains("errors"))  {
      drawOpt.ReplaceAll("errors","");
      DrawErrors(drawOpt);
   }
   else if (plotOpt.Contains("confidenceinterval") ||
            plotOpt.Contains("confinterval")) {
      // parse level option
      drawOpt.ReplaceAll("confidenceinterval","");
      drawOpt.ReplaceAll("confinterval","");
      Double_t level = 0.95;
      const char * s = strstr(plotOpt.Data(),"interval@");
      // coverity [secure_coding : FALSE]
      if (s != 0) sscanf(s,"interval@%lf",&level);
      if((level <= 0) || (level >= 1)) {
         Warning("Draw","given confidence level %.3lf is invalid - use default 0.95",level);
         level = 0.95;
      }
      DrawConfidenceInterval(drawOpt,level);
   }
   else {
      if (fPDF) delete fPDF;
      fPDF = GetKDEFunction();
      fPDF->Draw(drawOpt);
   }
}

void TKDE::DrawErrors(TString& drawOpt) {
   // Draws a TGraphErrors for the KDE errors
   if (fGraph) delete fGraph;
   fGraph = GetGraphWithErrors();
   fGraph->Draw(drawOpt.Data());
}

TGraphErrors* TKDE::GetGraphWithErrors(UInt_t npx, double xmin, double xmax) {
   if (xmin>= xmax) { xmin = fXMin; xmax = fXMax; }
   // return a TGraphErrors for the KDE errors
   UInt_t n = npx;
   Double_t* x = new Double_t[n + 1];
   Double_t* ex = new Double_t[n + 1];
   Double_t* y = new Double_t[n + 1];
   Double_t* ey = new Double_t[n + 1];
   for (UInt_t i = 0; i <= n; ++i) {
      x[i] = xmin + i * (xmax - xmin) / n;
      y[i] = (*this)(x[i]);
      ex[i] = 0;
      ey[i] = this->GetError(x[i]);
   }
   TGraphErrors* ge = new TGraphErrors(n, &x[0], &y[0], &ex[0], &ey[0]);
   ge->SetName("kde_graph_error");
   ge->SetTitle("Errors");
   delete [] x;
   delete [] ex;
   delete [] y;
   delete [] ey;
   return ge;
}

void TKDE::DrawConfidenceInterval(TString& drawOpt,double cl) {
   // Draws the KDE and its confidence interval
   GetKDEFunction()->Draw(drawOpt.Data());
   TF1* upper = GetPDFUpperConfidenceInterval(cl);
   upper->SetLineColor(kBlue);
   upper->Draw(("same" + drawOpt).Data());
   TF1* lower = GetPDFLowerConfidenceInterval(cl);
   lower->SetLineColor(kRed);
   lower->Draw(("same" + drawOpt).Data());
   if (fUpperPDF) delete fUpperPDF;
   if (fLowerPDF) delete fLowerPDF;
   fUpperPDF = upper;
   fLowerPDF = lower;
}

Double_t TKDE::GetFixedWeight() const {
   // Returns the bandwidth for the non adaptive KDE
   Double_t result = -1.;
   if (fIteration == TKDE::kAdaptive) {
      this->Warning("GetFixedWeight()", "Fixed iteration option not enabled. Returning %f.", result);
   } else {
      result = fKernel->GetFixedWeight();
   }
   return result;
}

const Double_t *  TKDE::GetAdaptiveWeights() const {
   // Returns the bandwidths for the adaptive KDE
   if (fIteration != TKDE::kAdaptive) {
      this->Warning("GetFixedWeight()", "Adaptive iteration option not enabled. Returning a NULL pointer<");
      return 0;
   }
   if (fNewData) (const_cast<TKDE*>(this))->InitFromNewData();
   return &(fKernel->GetAdaptiveWeights()).front();
}

Double_t TKDE::TKernel::GetFixedWeight() const {
   // Returns the bandwidth for the non adaptive KDE
   return fWeights[0];
}

const std::vector<Double_t> & TKDE::TKernel::GetAdaptiveWeights() const {
   // Returns the bandwidth for the non adaptive KDE
   return fWeights;
}

Double_t TKDE::TKernel::operator()(Double_t x) const {
   // The internal class's unary function: returns the kernel density estimate
   Double_t result(0.0);
   UInt_t n = fKDE->fData.size();
   // case of bins or weighted data 
   Bool_t useBins = (fKDE->fBinCount.size() == n);
   Double_t nSum = (useBins) ? fKDE->fSumOfCounts : fKDE->fNEvents;
   // double dmin = 1.E10;
   // double xmin,bmin,wmin; 
   for (UInt_t i = 0; i < n; ++i) {
      Double_t binCount = (useBins) ? fKDE->fBinCount[i] : 1.0;
      //printf("data point %i  %f  count %f weight % f result % f\n",i,fKDE->fData[i],binCount,fWeights[i], result);
      result += binCount / fWeights[i] * (*fKDE->fKernelFunction)((x - fKDE->fData[i]) / fWeights[i]);
      if (fKDE->fAsymLeft) {
         result -= binCount / fWeights[i] * (*fKDE->fKernelFunction)((x - (2. * fKDE->fXMin - fKDE->fData[i])) / fWeights[i]);
      }
      if (fKDE->fAsymRight) {
         result -= binCount / fWeights[i] * (*fKDE->fKernelFunction)((x - (2. * fKDE->fXMax - fKDE->fData[i])) / fWeights[i]);
      }
      // if ( TMath::IsNaN(result) ) {
      //    printf("event %i count %f  weight %f  data % f x %f \n",i,binCount,fWeights[i],fKDE->fData[i],x );
      // }
      // if ( result <= 0 ) {
      //    printf("event %i count %f  weight %f  data % f x %f \n",i,binCount,fWeights[i],fKDE->fData[i],x );
      // }
      // if (std::abs(x -  fKDE->fData[i]) < dmin ) {
      //    xmin = x;
      //    bmin = binCount;
      //    wmin = fWeights[i];
      //    dmin = std::abs(x -  fKDE->fData[i]);
      // }
   }
   if ( TMath::IsNaN(result) ) {
      fKDE->Warning("operator()","Result is NaN for  x %f \n",x);
    //xmin % f , %f, %f \n",result,x,xmin,bmin,wmin );
   }
   return result / nSum;
}

UInt_t TKDE::Index(Double_t x) const {
   // Returns the indices (bins) for the binned weights
   Int_t bin = Int_t((x - fXMin) * fWeightSize);
   if (bin == (Int_t)fData.size()) return --bin;
   if (fUseMirroring && (fMirrorLeft || !fMirrorRight)) {
      bin += fData.size() / (fMirrorLeft + fMirrorRight + 1);
   }
   if (bin > (Int_t)fData.size()) {
      return (Int_t)(fData.size()) - 1;
   } else if (bin <= 0) {
      return 0;
   }
   return bin;
}

Double_t TKDE::UpperConfidenceInterval(const Double_t* x, const Double_t* p) const {
   // Returns the pointwise upper estimated density
   Double_t f = (*this)(x);
   Double_t sigma = GetError(*x);
   Double_t prob = 1. - (1.-*p)/2;   // this is 1.-alpha/2
   Double_t z = ROOT::Math::normal_quantile(prob, 1.0);
   return f + z * sigma;
}

Double_t TKDE::LowerConfidenceInterval(const Double_t* x, const Double_t* p) const {
   // Returns the pointwise lower estimated density
   Double_t f = (*this)(x);
   Double_t sigma = GetError(*x);
   Double_t prob = (1.-*p)/2;    // this is alpha/2
   Double_t z = ROOT::Math::normal_quantile(prob, 1.0);
   return f + z * sigma;
}


Double_t TKDE::GetBias(Double_t x) const {
   // Returns the pointwise approximate estimated density bias
   ROOT::Math::WrappedFunction<const TKDE&> kern(*this);
   ROOT::Math::RichardsonDerivator rd;
   rd.SetFunction(kern);
   Double_t df2 = rd.Derivative2(x);
   Double_t weight = fKernel->GetWeight(x); // Bandwidth
   return  0.5 * fKernelSigmas2[fKernelType] * std::pow(weight, 2) * df2;
}
Double_t TKDE::GetError(Double_t x) const {
   // Returns the pointwise sigma of estimated density
   Double_t kernelL2Norm = ComputeKernelL2Norm();
   Double_t f = (*this)(x);
   Double_t weight = fKernel->GetWeight(x); // Bandwidth
   Double_t result = f * kernelL2Norm / (fNEvents * weight);
   return std::sqrt(result);
}

void TKDE::CheckKernelValidity() {
   // Checks if kernel has unit integral, mu = 0 and positive finite sigma conditions
   Double_t valid = kTRUE;
   Double_t unity = ComputeKernelIntegral();
   valid = valid && unity == 1.;
   if (!valid) {
      Error("CheckKernelValidity", "Kernel's integral is %f",unity);
   }
   Double_t mu = ComputeKernelMu();
   valid = valid && mu == 0.;
   if (!valid) {
      Error("CheckKernelValidity", "Kernel's mu is %f" ,mu);
   }
   Double_t sigma2 = ComputeKernelSigma2();
   valid = valid && sigma2 > 0 && sigma2 != std::numeric_limits<Double_t>::infinity();
   if (!valid) {
      Error("CheckKernelValidity", "Kernel's sigma2 is %f",sigma2);
   }
   if (!valid) {
      Error("CheckKernelValidity", "Validation conditions: the kernel's integral must be 1, the kernel's mu must be zero and the kernel's sigma2 must be finite positive to be a suitable kernel.");
      //exit(EXIT_FAILURE);
   }
}

Double_t TKDE::ComputeKernelL2Norm() const {
   // Computes the kernel's L2 norm
   ROOT::Math::IntegratorOneDim ig(ROOT::Math::IntegrationOneDim::kGAUSS);
   KernelIntegrand kernel(this, TKDE::KernelIntegrand::kNorm);
   ig.SetFunction(kernel);
   Double_t result = ig.Integral();
   return result;
}

Double_t TKDE::ComputeKernelSigma2() const {
   // Computes the kernel's sigma squared
   ROOT::Math::IntegratorOneDim ig( ROOT::Math::IntegrationOneDim::kGAUSS);
   KernelIntegrand kernel(this, TKDE::KernelIntegrand::kSigma2);
   ig.SetFunction(kernel);
   Double_t result = ig.Integral();
   return result;
}

Double_t TKDE::ComputeKernelMu() const {
   // Computes the kernel's mu
   ROOT::Math::IntegratorOneDim ig(ROOT::Math::IntegrationOneDim::kGAUSS);
   KernelIntegrand kernel(this, TKDE::KernelIntegrand::kMu);
   ig.SetFunction(kernel);
   Double_t result = ig.Integral();
   return result;
}

Double_t TKDE::ComputeKernelIntegral() const {
   // Computes the kernel's integral which ought to be unity
   ROOT::Math::IntegratorOneDim ig(ROOT::Math::IntegrationOneDim::kGAUSS);
   KernelIntegrand kernel(this, TKDE::KernelIntegrand::kUnitIntegration);
   ig.SetFunction(kernel);
   Double_t result = ig.Integral();
   return result;
}

void TKDE::ComputeDataStats() {
   /// in case of weights use
   if (!fEventWeights.empty() ) {
      // weighted data
      double x1 = fXMin - 0.001*(fXMax-fXMin);
      double x2 = fXMax + 0.001*(fXMax-fXMin);
      TH1D h1("temphist","", 500, x1, x2);
      h1.FillN(fEvents.size(), fEvents.data(), fEventWeights.data() );
      assert (h1.GetSumOfWeights() > 0) ;
      fMean = h1.GetMean();
      fSigma = h1.GetRMS();
      // compute robust sigma using midspread
      Double_t quantiles[2] = {0.0, 0.0};
      Double_t prob[2] = {0.25, 0.75};
      h1.GetQuantiles(2, quantiles, prob);
      Double_t midspread = quantiles[1] - quantiles[0];
      fSigmaRob = std::min(fSigma, midspread / 1.349); // Sigma's robust estimator
      //printf("weight case - stat: m = %f, s= %f, sr = %f \n",fMean, fSigma, midspread);
      return;
   }
   else {
      // compute statistics using the data
      SetMean();
      Double_t midspread = ComputeMidspread();
      SetSigma(midspread);
      //printf("un-weight case - stat: m = %f, s= %f, sr = %f \n",fMean, fSigma, fSigmaRob);
   }
}

Double_t TKDE::ComputeMidspread () {
   // Computes the inter-quartile range from the data
   std::sort(fEvents.begin(), fEvents.end());
   Double_t quantiles[2] = {0.0, 0.0};
   Double_t prob[2] = {0.25, 0.75};
   TMath::Quantiles(fEvents.size(), 2, &fEvents[0], quantiles, prob);
   Double_t lowerquartile = quantiles[0];
   Double_t upperquartile = quantiles[1];
   return upperquartile - lowerquartile;
}

void TKDE::SetUserCanonicalBandwidth() {
   // Computes the user's input kernel function canonical bandwidth
   fCanonicalBandwidths[kUserDefined] = std::pow(ComputeKernelL2Norm() / std::pow(ComputeKernelSigma2(), 2), 1. / 5.);
}

void TKDE::SetUserKernelSigma2() {
   // Computes the user's input kernel function sigma2
   fKernelSigmas2[kUserDefined] = ComputeKernelSigma2();
}

TKDE::KernelIntegrand::KernelIntegrand(const TKDE* kde, EIntegralResult intRes) : fKDE(kde), fIntegralResult(intRes) {
   // Internal class constructor
}

Double_t TKDE::KernelIntegrand::operator()(Double_t x) const {
   // Internal class unary function
   if (fIntegralResult == kNorm) {
      return std::pow((*fKDE->fKernelFunction)(x), 2);
   } else if (fIntegralResult == kMu) {
      return x * (*fKDE->fKernelFunction)(x);
   } else if (fIntegralResult == kSigma2) {
      return std::pow(x, 2) * (*fKDE->fKernelFunction)(x);
   } else if (fIntegralResult == kUnitIntegration) {
      return (*fKDE->fKernelFunction)(x);
   } else {
      return -1;
   }
}

TF1* TKDE::GetKDEFunction(UInt_t npx, Double_t xMin , Double_t xMax) {
   //Returns the estimated density
   TString name = "KDEFunc_"; name+= GetName();
   TString title = "KDE "; title+= GetTitle();
   if (xMin >= xMax) { xMin = fXMin; xMax = fXMax; }
   TF1 * pdf = new TF1(name.Data(), this, xMin, xMax, 0);
   if (npx > 0) pdf->SetNpx(npx);
   pdf->SetTitle(title);
   TF1 * f =  (TF1*)pdf->Clone();
   delete pdf;
   return f;
}

TF1* TKDE::GetPDFUpperConfidenceInterval(Double_t confidenceLevel, UInt_t npx, Double_t xMin , Double_t xMax) {
   // Returns the upper estimated density
   TString name;
   name.Form("KDE_UpperCL%f5.3_%s",confidenceLevel,GetName());
   if (xMin >= xMax) { xMin = fXMin; xMax = fXMax; }
   TF1 * upperPDF = new TF1(name, this, &TKDE::UpperConfidenceInterval, xMin, xMax, 1);
   upperPDF->SetParameter(0, confidenceLevel);
   if (npx > 0) upperPDF->SetNpx(npx);
   TF1 * f =  (TF1*)upperPDF->Clone();
   delete upperPDF;
   return f;
}

TF1* TKDE::GetPDFLowerConfidenceInterval(Double_t confidenceLevel, UInt_t npx, Double_t xMin , Double_t xMax) {
   // Returns the upper estimated density
   TString name;
   name.Form("KDE_LowerCL%f5.3_%s",confidenceLevel,GetName());
   if (xMin >= xMax) { xMin = fXMin; xMax = fXMax; }
   TF1 * lowerPDF = new TF1(name, this, &TKDE::LowerConfidenceInterval, xMin, xMax, 1);
   lowerPDF->SetParameter(0, confidenceLevel);
   if (npx > 0) lowerPDF->SetNpx(npx);
   TF1 * f = (TF1*)lowerPDF->Clone();
   delete lowerPDF;
   return f;
}

TF1* TKDE::GetKDEApproximateBias(UInt_t npx, Double_t xMin , Double_t xMax) {
   // Returns the approximate bias
   TString name = "KDE_Bias_"; name += GetName();
   if (xMin >= xMax) { xMin = fXMin; xMax = fXMax; }
   TF1 * approximateBias = new TF1(name, this, &TKDE::ApproximateBias, xMin, xMax, 0);
   if (npx > 0) approximateBias->SetNpx(npx);
   TF1 * f =  (TF1*)approximateBias->Clone();
   delete approximateBias;
   return f;
}
 TKDE.cxx:1
 TKDE.cxx:2
 TKDE.cxx:3
 TKDE.cxx:4
 TKDE.cxx:5
 TKDE.cxx:6
 TKDE.cxx:7
 TKDE.cxx:8
 TKDE.cxx:9
 TKDE.cxx:10
 TKDE.cxx:11
 TKDE.cxx:12
 TKDE.cxx:13
 TKDE.cxx:14
 TKDE.cxx:15
 TKDE.cxx:16
 TKDE.cxx:17
 TKDE.cxx:18
 TKDE.cxx:19
 TKDE.cxx:20
 TKDE.cxx:21
 TKDE.cxx:22
 TKDE.cxx:23
 TKDE.cxx:24
 TKDE.cxx:25
 TKDE.cxx:26
 TKDE.cxx:27
 TKDE.cxx:28
 TKDE.cxx:29
 TKDE.cxx:30
 TKDE.cxx:31
 TKDE.cxx:32
 TKDE.cxx:33
 TKDE.cxx:34
 TKDE.cxx:35
 TKDE.cxx:36
 TKDE.cxx:37
 TKDE.cxx:38
 TKDE.cxx:39
 TKDE.cxx:40
 TKDE.cxx:41
 TKDE.cxx:42
 TKDE.cxx:43
 TKDE.cxx:44
 TKDE.cxx:45
 TKDE.cxx:46
 TKDE.cxx:47
 TKDE.cxx:48
 TKDE.cxx:49
 TKDE.cxx:50
 TKDE.cxx:51
 TKDE.cxx:52
 TKDE.cxx:53
 TKDE.cxx:54
 TKDE.cxx:55
 TKDE.cxx:56
 TKDE.cxx:57
 TKDE.cxx:58
 TKDE.cxx:59
 TKDE.cxx:60
 TKDE.cxx:61
 TKDE.cxx:62
 TKDE.cxx:63
 TKDE.cxx:64
 TKDE.cxx:65
 TKDE.cxx:66
 TKDE.cxx:67
 TKDE.cxx:68
 TKDE.cxx:69
 TKDE.cxx:70
 TKDE.cxx:71
 TKDE.cxx:72
 TKDE.cxx:73
 TKDE.cxx:74
 TKDE.cxx:75
 TKDE.cxx:76
 TKDE.cxx:77
 TKDE.cxx:78
 TKDE.cxx:79
 TKDE.cxx:80
 TKDE.cxx:81
 TKDE.cxx:82
 TKDE.cxx:83
 TKDE.cxx:84
 TKDE.cxx:85
 TKDE.cxx:86
 TKDE.cxx:87
 TKDE.cxx:88
 TKDE.cxx:89
 TKDE.cxx:90
 TKDE.cxx:91
 TKDE.cxx:92
 TKDE.cxx:93
 TKDE.cxx:94
 TKDE.cxx:95
 TKDE.cxx:96
 TKDE.cxx:97
 TKDE.cxx:98
 TKDE.cxx:99
 TKDE.cxx:100
 TKDE.cxx:101
 TKDE.cxx:102
 TKDE.cxx:103
 TKDE.cxx:104
 TKDE.cxx:105
 TKDE.cxx:106
 TKDE.cxx:107
 TKDE.cxx:108
 TKDE.cxx:109
 TKDE.cxx:110
 TKDE.cxx:111
 TKDE.cxx:112
 TKDE.cxx:113
 TKDE.cxx:114
 TKDE.cxx:115
 TKDE.cxx:116
 TKDE.cxx:117
 TKDE.cxx:118
 TKDE.cxx:119
 TKDE.cxx:120
 TKDE.cxx:121
 TKDE.cxx:122
 TKDE.cxx:123
 TKDE.cxx:124
 TKDE.cxx:125
 TKDE.cxx:126
 TKDE.cxx:127
 TKDE.cxx:128
 TKDE.cxx:129
 TKDE.cxx:130
 TKDE.cxx:131
 TKDE.cxx:132
 TKDE.cxx:133
 TKDE.cxx:134
 TKDE.cxx:135
 TKDE.cxx:136
 TKDE.cxx:137
 TKDE.cxx:138
 TKDE.cxx:139
 TKDE.cxx:140
 TKDE.cxx:141
 TKDE.cxx:142
 TKDE.cxx:143
 TKDE.cxx:144
 TKDE.cxx:145
 TKDE.cxx:146
 TKDE.cxx:147
 TKDE.cxx:148
 TKDE.cxx:149
 TKDE.cxx:150
 TKDE.cxx:151
 TKDE.cxx:152
 TKDE.cxx:153
 TKDE.cxx:154
 TKDE.cxx:155
 TKDE.cxx:156
 TKDE.cxx:157
 TKDE.cxx:158
 TKDE.cxx:159
 TKDE.cxx:160
 TKDE.cxx:161
 TKDE.cxx:162
 TKDE.cxx:163
 TKDE.cxx:164
 TKDE.cxx:165
 TKDE.cxx:166
 TKDE.cxx:167
 TKDE.cxx:168
 TKDE.cxx:169
 TKDE.cxx:170
 TKDE.cxx:171
 TKDE.cxx:172
 TKDE.cxx:173
 TKDE.cxx:174
 TKDE.cxx:175
 TKDE.cxx:176
 TKDE.cxx:177
 TKDE.cxx:178
 TKDE.cxx:179
 TKDE.cxx:180
 TKDE.cxx:181
 TKDE.cxx:182
 TKDE.cxx:183
 TKDE.cxx:184
 TKDE.cxx:185
 TKDE.cxx:186
 TKDE.cxx:187
 TKDE.cxx:188
 TKDE.cxx:189
 TKDE.cxx:190
 TKDE.cxx:191
 TKDE.cxx:192
 TKDE.cxx:193
 TKDE.cxx:194
 TKDE.cxx:195
 TKDE.cxx:196
 TKDE.cxx:197
 TKDE.cxx:198
 TKDE.cxx:199
 TKDE.cxx:200
 TKDE.cxx:201
 TKDE.cxx:202
 TKDE.cxx:203
 TKDE.cxx:204
 TKDE.cxx:205
 TKDE.cxx:206
 TKDE.cxx:207
 TKDE.cxx:208
 TKDE.cxx:209
 TKDE.cxx:210
 TKDE.cxx:211
 TKDE.cxx:212
 TKDE.cxx:213
 TKDE.cxx:214
 TKDE.cxx:215
 TKDE.cxx:216
 TKDE.cxx:217
 TKDE.cxx:218
 TKDE.cxx:219
 TKDE.cxx:220
 TKDE.cxx:221
 TKDE.cxx:222
 TKDE.cxx:223
 TKDE.cxx:224
 TKDE.cxx:225
 TKDE.cxx:226
 TKDE.cxx:227
 TKDE.cxx:228
 TKDE.cxx:229
 TKDE.cxx:230
 TKDE.cxx:231
 TKDE.cxx:232
 TKDE.cxx:233
 TKDE.cxx:234
 TKDE.cxx:235
 TKDE.cxx:236
 TKDE.cxx:237
 TKDE.cxx:238
 TKDE.cxx:239
 TKDE.cxx:240
 TKDE.cxx:241
 TKDE.cxx:242
 TKDE.cxx:243
 TKDE.cxx:244
 TKDE.cxx:245
 TKDE.cxx:246
 TKDE.cxx:247
 TKDE.cxx:248
 TKDE.cxx:249
 TKDE.cxx:250
 TKDE.cxx:251
 TKDE.cxx:252
 TKDE.cxx:253
 TKDE.cxx:254
 TKDE.cxx:255
 TKDE.cxx:256
 TKDE.cxx:257
 TKDE.cxx:258
 TKDE.cxx:259
 TKDE.cxx:260
 TKDE.cxx:261
 TKDE.cxx:262
 TKDE.cxx:263
 TKDE.cxx:264
 TKDE.cxx:265
 TKDE.cxx:266
 TKDE.cxx:267
 TKDE.cxx:268
 TKDE.cxx:269
 TKDE.cxx:270
 TKDE.cxx:271
 TKDE.cxx:272
 TKDE.cxx:273
 TKDE.cxx:274
 TKDE.cxx:275
 TKDE.cxx:276
 TKDE.cxx:277
 TKDE.cxx:278
 TKDE.cxx:279
 TKDE.cxx:280
 TKDE.cxx:281
 TKDE.cxx:282
 TKDE.cxx:283
 TKDE.cxx:284
 TKDE.cxx:285
 TKDE.cxx:286
 TKDE.cxx:287
 TKDE.cxx:288
 TKDE.cxx:289
 TKDE.cxx:290
 TKDE.cxx:291
 TKDE.cxx:292
 TKDE.cxx:293
 TKDE.cxx:294
 TKDE.cxx:295
 TKDE.cxx:296
 TKDE.cxx:297
 TKDE.cxx:298
 TKDE.cxx:299
 TKDE.cxx:300
 TKDE.cxx:301
 TKDE.cxx:302
 TKDE.cxx:303
 TKDE.cxx:304
 TKDE.cxx:305
 TKDE.cxx:306
 TKDE.cxx:307
 TKDE.cxx:308
 TKDE.cxx:309
 TKDE.cxx:310
 TKDE.cxx:311
 TKDE.cxx:312
 TKDE.cxx:313
 TKDE.cxx:314
 TKDE.cxx:315
 TKDE.cxx:316
 TKDE.cxx:317
 TKDE.cxx:318
 TKDE.cxx:319
 TKDE.cxx:320
 TKDE.cxx:321
 TKDE.cxx:322
 TKDE.cxx:323
 TKDE.cxx:324
 TKDE.cxx:325
 TKDE.cxx:326
 TKDE.cxx:327
 TKDE.cxx:328
 TKDE.cxx:329
 TKDE.cxx:330
 TKDE.cxx:331
 TKDE.cxx:332
 TKDE.cxx:333
 TKDE.cxx:334
 TKDE.cxx:335
 TKDE.cxx:336
 TKDE.cxx:337
 TKDE.cxx:338
 TKDE.cxx:339
 TKDE.cxx:340
 TKDE.cxx:341
 TKDE.cxx:342
 TKDE.cxx:343
 TKDE.cxx:344
 TKDE.cxx:345
 TKDE.cxx:346
 TKDE.cxx:347
 TKDE.cxx:348
 TKDE.cxx:349
 TKDE.cxx:350
 TKDE.cxx:351
 TKDE.cxx:352
 TKDE.cxx:353
 TKDE.cxx:354
 TKDE.cxx:355
 TKDE.cxx:356
 TKDE.cxx:357
 TKDE.cxx:358
 TKDE.cxx:359
 TKDE.cxx:360
 TKDE.cxx:361
 TKDE.cxx:362
 TKDE.cxx:363
 TKDE.cxx:364
 TKDE.cxx:365
 TKDE.cxx:366
 TKDE.cxx:367
 TKDE.cxx:368
 TKDE.cxx:369
 TKDE.cxx:370
 TKDE.cxx:371
 TKDE.cxx:372
 TKDE.cxx:373
 TKDE.cxx:374
 TKDE.cxx:375
 TKDE.cxx:376
 TKDE.cxx:377
 TKDE.cxx:378
 TKDE.cxx:379
 TKDE.cxx:380
 TKDE.cxx:381
 TKDE.cxx:382
 TKDE.cxx:383
 TKDE.cxx:384
 TKDE.cxx:385
 TKDE.cxx:386
 TKDE.cxx:387
 TKDE.cxx:388
 TKDE.cxx:389
 TKDE.cxx:390
 TKDE.cxx:391
 TKDE.cxx:392
 TKDE.cxx:393
 TKDE.cxx:394
 TKDE.cxx:395
 TKDE.cxx:396
 TKDE.cxx:397
 TKDE.cxx:398
 TKDE.cxx:399
 TKDE.cxx:400
 TKDE.cxx:401
 TKDE.cxx:402
 TKDE.cxx:403
 TKDE.cxx:404
 TKDE.cxx:405
 TKDE.cxx:406
 TKDE.cxx:407
 TKDE.cxx:408
 TKDE.cxx:409
 TKDE.cxx:410
 TKDE.cxx:411
 TKDE.cxx:412
 TKDE.cxx:413
 TKDE.cxx:414
 TKDE.cxx:415
 TKDE.cxx:416
 TKDE.cxx:417
 TKDE.cxx:418
 TKDE.cxx:419
 TKDE.cxx:420
 TKDE.cxx:421
 TKDE.cxx:422
 TKDE.cxx:423
 TKDE.cxx:424
 TKDE.cxx:425
 TKDE.cxx:426
 TKDE.cxx:427
 TKDE.cxx:428
 TKDE.cxx:429
 TKDE.cxx:430
 TKDE.cxx:431
 TKDE.cxx:432
 TKDE.cxx:433
 TKDE.cxx:434
 TKDE.cxx:435
 TKDE.cxx:436
 TKDE.cxx:437
 TKDE.cxx:438
 TKDE.cxx:439
 TKDE.cxx:440
 TKDE.cxx:441
 TKDE.cxx:442
 TKDE.cxx:443
 TKDE.cxx:444
 TKDE.cxx:445
 TKDE.cxx:446
 TKDE.cxx:447
 TKDE.cxx:448
 TKDE.cxx:449
 TKDE.cxx:450
 TKDE.cxx:451
 TKDE.cxx:452
 TKDE.cxx:453
 TKDE.cxx:454
 TKDE.cxx:455
 TKDE.cxx:456
 TKDE.cxx:457
 TKDE.cxx:458
 TKDE.cxx:459
 TKDE.cxx:460
 TKDE.cxx:461
 TKDE.cxx:462
 TKDE.cxx:463
 TKDE.cxx:464
 TKDE.cxx:465
 TKDE.cxx:466
 TKDE.cxx:467
 TKDE.cxx:468
 TKDE.cxx:469
 TKDE.cxx:470
 TKDE.cxx:471
 TKDE.cxx:472
 TKDE.cxx:473
 TKDE.cxx:474
 TKDE.cxx:475
 TKDE.cxx:476
 TKDE.cxx:477
 TKDE.cxx:478
 TKDE.cxx:479
 TKDE.cxx:480
 TKDE.cxx:481
 TKDE.cxx:482
 TKDE.cxx:483
 TKDE.cxx:484
 TKDE.cxx:485
 TKDE.cxx:486
 TKDE.cxx:487
 TKDE.cxx:488
 TKDE.cxx:489
 TKDE.cxx:490
 TKDE.cxx:491
 TKDE.cxx:492
 TKDE.cxx:493
 TKDE.cxx:494
 TKDE.cxx:495
 TKDE.cxx:496
 TKDE.cxx:497
 TKDE.cxx:498
 TKDE.cxx:499
 TKDE.cxx:500
 TKDE.cxx:501
 TKDE.cxx:502
 TKDE.cxx:503
 TKDE.cxx:504
 TKDE.cxx:505
 TKDE.cxx:506
 TKDE.cxx:507
 TKDE.cxx:508
 TKDE.cxx:509
 TKDE.cxx:510
 TKDE.cxx:511
 TKDE.cxx:512
 TKDE.cxx:513
 TKDE.cxx:514
 TKDE.cxx:515
 TKDE.cxx:516
 TKDE.cxx:517
 TKDE.cxx:518
 TKDE.cxx:519
 TKDE.cxx:520
 TKDE.cxx:521
 TKDE.cxx:522
 TKDE.cxx:523
 TKDE.cxx:524
 TKDE.cxx:525
 TKDE.cxx:526
 TKDE.cxx:527
 TKDE.cxx:528
 TKDE.cxx:529
 TKDE.cxx:530
 TKDE.cxx:531
 TKDE.cxx:532
 TKDE.cxx:533
 TKDE.cxx:534
 TKDE.cxx:535
 TKDE.cxx:536
 TKDE.cxx:537
 TKDE.cxx:538
 TKDE.cxx:539
 TKDE.cxx:540
 TKDE.cxx:541
 TKDE.cxx:542
 TKDE.cxx:543
 TKDE.cxx:544
 TKDE.cxx:545
 TKDE.cxx:546
 TKDE.cxx:547
 TKDE.cxx:548
 TKDE.cxx:549
 TKDE.cxx:550
 TKDE.cxx:551
 TKDE.cxx:552
 TKDE.cxx:553
 TKDE.cxx:554
 TKDE.cxx:555
 TKDE.cxx:556
 TKDE.cxx:557
 TKDE.cxx:558
 TKDE.cxx:559
 TKDE.cxx:560
 TKDE.cxx:561
 TKDE.cxx:562
 TKDE.cxx:563
 TKDE.cxx:564
 TKDE.cxx:565
 TKDE.cxx:566
 TKDE.cxx:567
 TKDE.cxx:568
 TKDE.cxx:569
 TKDE.cxx:570
 TKDE.cxx:571
 TKDE.cxx:572
 TKDE.cxx:573
 TKDE.cxx:574
 TKDE.cxx:575
 TKDE.cxx:576
 TKDE.cxx:577
 TKDE.cxx:578
 TKDE.cxx:579
 TKDE.cxx:580
 TKDE.cxx:581
 TKDE.cxx:582
 TKDE.cxx:583
 TKDE.cxx:584
 TKDE.cxx:585
 TKDE.cxx:586
 TKDE.cxx:587
 TKDE.cxx:588
 TKDE.cxx:589
 TKDE.cxx:590
 TKDE.cxx:591
 TKDE.cxx:592
 TKDE.cxx:593
 TKDE.cxx:594
 TKDE.cxx:595
 TKDE.cxx:596
 TKDE.cxx:597
 TKDE.cxx:598
 TKDE.cxx:599
 TKDE.cxx:600
 TKDE.cxx:601
 TKDE.cxx:602
 TKDE.cxx:603
 TKDE.cxx:604
 TKDE.cxx:605
 TKDE.cxx:606
 TKDE.cxx:607
 TKDE.cxx:608
 TKDE.cxx:609
 TKDE.cxx:610
 TKDE.cxx:611
 TKDE.cxx:612
 TKDE.cxx:613
 TKDE.cxx:614
 TKDE.cxx:615
 TKDE.cxx:616
 TKDE.cxx:617
 TKDE.cxx:618
 TKDE.cxx:619
 TKDE.cxx:620
 TKDE.cxx:621
 TKDE.cxx:622
 TKDE.cxx:623
 TKDE.cxx:624
 TKDE.cxx:625
 TKDE.cxx:626
 TKDE.cxx:627
 TKDE.cxx:628
 TKDE.cxx:629
 TKDE.cxx:630
 TKDE.cxx:631
 TKDE.cxx:632
 TKDE.cxx:633
 TKDE.cxx:634
 TKDE.cxx:635
 TKDE.cxx:636
 TKDE.cxx:637
 TKDE.cxx:638
 TKDE.cxx:639
 TKDE.cxx:640
 TKDE.cxx:641
 TKDE.cxx:642
 TKDE.cxx:643
 TKDE.cxx:644
 TKDE.cxx:645
 TKDE.cxx:646
 TKDE.cxx:647
 TKDE.cxx:648
 TKDE.cxx:649
 TKDE.cxx:650
 TKDE.cxx:651
 TKDE.cxx:652
 TKDE.cxx:653
 TKDE.cxx:654
 TKDE.cxx:655
 TKDE.cxx:656
 TKDE.cxx:657
 TKDE.cxx:658
 TKDE.cxx:659
 TKDE.cxx:660
 TKDE.cxx:661
 TKDE.cxx:662
 TKDE.cxx:663
 TKDE.cxx:664
 TKDE.cxx:665
 TKDE.cxx:666
 TKDE.cxx:667
 TKDE.cxx:668
 TKDE.cxx:669
 TKDE.cxx:670
 TKDE.cxx:671
 TKDE.cxx:672
 TKDE.cxx:673
 TKDE.cxx:674
 TKDE.cxx:675
 TKDE.cxx:676
 TKDE.cxx:677
 TKDE.cxx:678
 TKDE.cxx:679
 TKDE.cxx:680
 TKDE.cxx:681
 TKDE.cxx:682
 TKDE.cxx:683
 TKDE.cxx:684
 TKDE.cxx:685
 TKDE.cxx:686
 TKDE.cxx:687
 TKDE.cxx:688
 TKDE.cxx:689
 TKDE.cxx:690
 TKDE.cxx:691
 TKDE.cxx:692
 TKDE.cxx:693
 TKDE.cxx:694
 TKDE.cxx:695
 TKDE.cxx:696
 TKDE.cxx:697
 TKDE.cxx:698
 TKDE.cxx:699
 TKDE.cxx:700
 TKDE.cxx:701
 TKDE.cxx:702
 TKDE.cxx:703
 TKDE.cxx:704
 TKDE.cxx:705
 TKDE.cxx:706
 TKDE.cxx:707
 TKDE.cxx:708
 TKDE.cxx:709
 TKDE.cxx:710
 TKDE.cxx:711
 TKDE.cxx:712
 TKDE.cxx:713
 TKDE.cxx:714
 TKDE.cxx:715
 TKDE.cxx:716
 TKDE.cxx:717
 TKDE.cxx:718
 TKDE.cxx:719
 TKDE.cxx:720
 TKDE.cxx:721
 TKDE.cxx:722
 TKDE.cxx:723
 TKDE.cxx:724
 TKDE.cxx:725
 TKDE.cxx:726
 TKDE.cxx:727
 TKDE.cxx:728
 TKDE.cxx:729
 TKDE.cxx:730
 TKDE.cxx:731
 TKDE.cxx:732
 TKDE.cxx:733
 TKDE.cxx:734
 TKDE.cxx:735
 TKDE.cxx:736
 TKDE.cxx:737
 TKDE.cxx:738
 TKDE.cxx:739
 TKDE.cxx:740
 TKDE.cxx:741
 TKDE.cxx:742
 TKDE.cxx:743
 TKDE.cxx:744
 TKDE.cxx:745
 TKDE.cxx:746
 TKDE.cxx:747
 TKDE.cxx:748
 TKDE.cxx:749
 TKDE.cxx:750
 TKDE.cxx:751
 TKDE.cxx:752
 TKDE.cxx:753
 TKDE.cxx:754
 TKDE.cxx:755
 TKDE.cxx:756
 TKDE.cxx:757
 TKDE.cxx:758
 TKDE.cxx:759
 TKDE.cxx:760
 TKDE.cxx:761
 TKDE.cxx:762
 TKDE.cxx:763
 TKDE.cxx:764
 TKDE.cxx:765
 TKDE.cxx:766
 TKDE.cxx:767
 TKDE.cxx:768
 TKDE.cxx:769
 TKDE.cxx:770
 TKDE.cxx:771
 TKDE.cxx:772
 TKDE.cxx:773
 TKDE.cxx:774
 TKDE.cxx:775
 TKDE.cxx:776
 TKDE.cxx:777
 TKDE.cxx:778
 TKDE.cxx:779
 TKDE.cxx:780
 TKDE.cxx:781
 TKDE.cxx:782
 TKDE.cxx:783
 TKDE.cxx:784
 TKDE.cxx:785
 TKDE.cxx:786
 TKDE.cxx:787
 TKDE.cxx:788
 TKDE.cxx:789
 TKDE.cxx:790
 TKDE.cxx:791
 TKDE.cxx:792
 TKDE.cxx:793
 TKDE.cxx:794
 TKDE.cxx:795
 TKDE.cxx:796
 TKDE.cxx:797
 TKDE.cxx:798
 TKDE.cxx:799
 TKDE.cxx:800
 TKDE.cxx:801
 TKDE.cxx:802
 TKDE.cxx:803
 TKDE.cxx:804
 TKDE.cxx:805
 TKDE.cxx:806
 TKDE.cxx:807
 TKDE.cxx:808
 TKDE.cxx:809
 TKDE.cxx:810
 TKDE.cxx:811
 TKDE.cxx:812
 TKDE.cxx:813
 TKDE.cxx:814
 TKDE.cxx:815
 TKDE.cxx:816
 TKDE.cxx:817
 TKDE.cxx:818
 TKDE.cxx:819
 TKDE.cxx:820
 TKDE.cxx:821
 TKDE.cxx:822
 TKDE.cxx:823
 TKDE.cxx:824
 TKDE.cxx:825
 TKDE.cxx:826
 TKDE.cxx:827
 TKDE.cxx:828
 TKDE.cxx:829
 TKDE.cxx:830
 TKDE.cxx:831
 TKDE.cxx:832
 TKDE.cxx:833
 TKDE.cxx:834
 TKDE.cxx:835
 TKDE.cxx:836
 TKDE.cxx:837
 TKDE.cxx:838
 TKDE.cxx:839
 TKDE.cxx:840
 TKDE.cxx:841
 TKDE.cxx:842
 TKDE.cxx:843
 TKDE.cxx:844
 TKDE.cxx:845
 TKDE.cxx:846
 TKDE.cxx:847
 TKDE.cxx:848
 TKDE.cxx:849
 TKDE.cxx:850
 TKDE.cxx:851
 TKDE.cxx:852
 TKDE.cxx:853
 TKDE.cxx:854
 TKDE.cxx:855
 TKDE.cxx:856
 TKDE.cxx:857
 TKDE.cxx:858
 TKDE.cxx:859
 TKDE.cxx:860
 TKDE.cxx:861
 TKDE.cxx:862
 TKDE.cxx:863
 TKDE.cxx:864
 TKDE.cxx:865
 TKDE.cxx:866
 TKDE.cxx:867
 TKDE.cxx:868
 TKDE.cxx:869
 TKDE.cxx:870
 TKDE.cxx:871
 TKDE.cxx:872
 TKDE.cxx:873
 TKDE.cxx:874
 TKDE.cxx:875
 TKDE.cxx:876
 TKDE.cxx:877
 TKDE.cxx:878
 TKDE.cxx:879
 TKDE.cxx:880
 TKDE.cxx:881
 TKDE.cxx:882
 TKDE.cxx:883
 TKDE.cxx:884
 TKDE.cxx:885
 TKDE.cxx:886
 TKDE.cxx:887
 TKDE.cxx:888
 TKDE.cxx:889
 TKDE.cxx:890
 TKDE.cxx:891
 TKDE.cxx:892
 TKDE.cxx:893
 TKDE.cxx:894
 TKDE.cxx:895
 TKDE.cxx:896
 TKDE.cxx:897
 TKDE.cxx:898
 TKDE.cxx:899
 TKDE.cxx:900
 TKDE.cxx:901
 TKDE.cxx:902
 TKDE.cxx:903
 TKDE.cxx:904
 TKDE.cxx:905
 TKDE.cxx:906
 TKDE.cxx:907
 TKDE.cxx:908
 TKDE.cxx:909
 TKDE.cxx:910
 TKDE.cxx:911
 TKDE.cxx:912
 TKDE.cxx:913
 TKDE.cxx:914
 TKDE.cxx:915
 TKDE.cxx:916
 TKDE.cxx:917
 TKDE.cxx:918
 TKDE.cxx:919
 TKDE.cxx:920
 TKDE.cxx:921
 TKDE.cxx:922
 TKDE.cxx:923
 TKDE.cxx:924
 TKDE.cxx:925
 TKDE.cxx:926
 TKDE.cxx:927
 TKDE.cxx:928
 TKDE.cxx:929
 TKDE.cxx:930
 TKDE.cxx:931
 TKDE.cxx:932
 TKDE.cxx:933
 TKDE.cxx:934
 TKDE.cxx:935
 TKDE.cxx:936
 TKDE.cxx:937
 TKDE.cxx:938
 TKDE.cxx:939
 TKDE.cxx:940
 TKDE.cxx:941
 TKDE.cxx:942
 TKDE.cxx:943
 TKDE.cxx:944
 TKDE.cxx:945
 TKDE.cxx:946
 TKDE.cxx:947
 TKDE.cxx:948
 TKDE.cxx:949
 TKDE.cxx:950
 TKDE.cxx:951
 TKDE.cxx:952
 TKDE.cxx:953
 TKDE.cxx:954
 TKDE.cxx:955
 TKDE.cxx:956
 TKDE.cxx:957
 TKDE.cxx:958
 TKDE.cxx:959
 TKDE.cxx:960
 TKDE.cxx:961
 TKDE.cxx:962
 TKDE.cxx:963
 TKDE.cxx:964
 TKDE.cxx:965
 TKDE.cxx:966
 TKDE.cxx:967
 TKDE.cxx:968
 TKDE.cxx:969
 TKDE.cxx:970
 TKDE.cxx:971
 TKDE.cxx:972
 TKDE.cxx:973
 TKDE.cxx:974
 TKDE.cxx:975
 TKDE.cxx:976
 TKDE.cxx:977
 TKDE.cxx:978
 TKDE.cxx:979
 TKDE.cxx:980
 TKDE.cxx:981
 TKDE.cxx:982
 TKDE.cxx:983
 TKDE.cxx:984
 TKDE.cxx:985
 TKDE.cxx:986
 TKDE.cxx:987
 TKDE.cxx:988
 TKDE.cxx:989
 TKDE.cxx:990
 TKDE.cxx:991
 TKDE.cxx:992
 TKDE.cxx:993
 TKDE.cxx:994
 TKDE.cxx:995
 TKDE.cxx:996
 TKDE.cxx:997
 TKDE.cxx:998
 TKDE.cxx:999
 TKDE.cxx:1000
 TKDE.cxx:1001
 TKDE.cxx:1002
 TKDE.cxx:1003
 TKDE.cxx:1004
 TKDE.cxx:1005
 TKDE.cxx:1006
 TKDE.cxx:1007
 TKDE.cxx:1008
 TKDE.cxx:1009
 TKDE.cxx:1010
 TKDE.cxx:1011
 TKDE.cxx:1012
 TKDE.cxx:1013
 TKDE.cxx:1014
 TKDE.cxx:1015
 TKDE.cxx:1016
 TKDE.cxx:1017
 TKDE.cxx:1018
 TKDE.cxx:1019
 TKDE.cxx:1020
 TKDE.cxx:1021
 TKDE.cxx:1022
 TKDE.cxx:1023
 TKDE.cxx:1024
 TKDE.cxx:1025
 TKDE.cxx:1026
 TKDE.cxx:1027
 TKDE.cxx:1028
 TKDE.cxx:1029
 TKDE.cxx:1030
 TKDE.cxx:1031
 TKDE.cxx:1032
 TKDE.cxx:1033
 TKDE.cxx:1034
 TKDE.cxx:1035
 TKDE.cxx:1036
 TKDE.cxx:1037
 TKDE.cxx:1038
 TKDE.cxx:1039
 TKDE.cxx:1040
 TKDE.cxx:1041
 TKDE.cxx:1042
 TKDE.cxx:1043
 TKDE.cxx:1044
 TKDE.cxx:1045
 TKDE.cxx:1046
 TKDE.cxx:1047
 TKDE.cxx:1048
 TKDE.cxx:1049
 TKDE.cxx:1050
 TKDE.cxx:1051
 TKDE.cxx:1052
 TKDE.cxx:1053
 TKDE.cxx:1054
 TKDE.cxx:1055
 TKDE.cxx:1056
 TKDE.cxx:1057
 TKDE.cxx:1058
 TKDE.cxx:1059
 TKDE.cxx:1060
 TKDE.cxx:1061
 TKDE.cxx:1062
 TKDE.cxx:1063
 TKDE.cxx:1064
 TKDE.cxx:1065
 TKDE.cxx:1066
 TKDE.cxx:1067
 TKDE.cxx:1068
 TKDE.cxx:1069
 TKDE.cxx:1070
 TKDE.cxx:1071
 TKDE.cxx:1072
 TKDE.cxx:1073
 TKDE.cxx:1074
 TKDE.cxx:1075
 TKDE.cxx:1076
 TKDE.cxx:1077
 TKDE.cxx:1078
 TKDE.cxx:1079
 TKDE.cxx:1080
 TKDE.cxx:1081
 TKDE.cxx:1082
 TKDE.cxx:1083
 TKDE.cxx:1084
 TKDE.cxx:1085
 TKDE.cxx:1086
 TKDE.cxx:1087
 TKDE.cxx:1088
 TKDE.cxx:1089
 TKDE.cxx:1090
 TKDE.cxx:1091
 TKDE.cxx:1092
 TKDE.cxx:1093
 TKDE.cxx:1094
 TKDE.cxx:1095
 TKDE.cxx:1096
 TKDE.cxx:1097
 TKDE.cxx:1098
 TKDE.cxx:1099
 TKDE.cxx:1100
 TKDE.cxx:1101
 TKDE.cxx:1102
 TKDE.cxx:1103
 TKDE.cxx:1104
 TKDE.cxx:1105
 TKDE.cxx:1106
 TKDE.cxx:1107
 TKDE.cxx:1108
 TKDE.cxx:1109
 TKDE.cxx:1110
 TKDE.cxx:1111
 TKDE.cxx:1112
 TKDE.cxx:1113
 TKDE.cxx:1114
 TKDE.cxx:1115
 TKDE.cxx:1116
 TKDE.cxx:1117
 TKDE.cxx:1118
 TKDE.cxx:1119
 TKDE.cxx:1120
 TKDE.cxx:1121
 TKDE.cxx:1122
 TKDE.cxx:1123
 TKDE.cxx:1124
 TKDE.cxx:1125
 TKDE.cxx:1126
 TKDE.cxx:1127
 TKDE.cxx:1128
 TKDE.cxx:1129
 TKDE.cxx:1130
 TKDE.cxx:1131
 TKDE.cxx:1132
 TKDE.cxx:1133
 TKDE.cxx:1134
 TKDE.cxx:1135
 TKDE.cxx:1136
 TKDE.cxx:1137
 TKDE.cxx:1138
 TKDE.cxx:1139
 TKDE.cxx:1140
 TKDE.cxx:1141
 TKDE.cxx:1142
 TKDE.cxx:1143
 TKDE.cxx:1144
 TKDE.cxx:1145
 TKDE.cxx:1146
 TKDE.cxx:1147
 TKDE.cxx:1148
 TKDE.cxx:1149
 TKDE.cxx:1150
 TKDE.cxx:1151
 TKDE.cxx:1152
 TKDE.cxx:1153
 TKDE.cxx:1154
 TKDE.cxx:1155
 TKDE.cxx:1156
 TKDE.cxx:1157
 TKDE.cxx:1158
 TKDE.cxx:1159
 TKDE.cxx:1160
 TKDE.cxx:1161
 TKDE.cxx:1162
 TKDE.cxx:1163
 TKDE.cxx:1164
 TKDE.cxx:1165
 TKDE.cxx:1166
 TKDE.cxx:1167
 TKDE.cxx:1168
 TKDE.cxx:1169