ROOT logo
// @(#)root/tmva $Id: MethodCFMlpANN.cxx 29195 2009-06-24 10:39:49Z brun $    
// Author: Andreas Hoecker, Joerg Stelzer, Helge Voss, Kai Voss 

/**********************************************************************************
 * Project: TMVA - a Root-integrated toolkit for multivariate Data analysis       *
 * Package: TMVA                                                                  *
 * Class  : TMVA::MethodCFMlpANN                                                  *
 * Web    : http://tmva.sourceforge.net                                           *
 *                                                                                *
 * Description:                                                                   *
 *      Implementation (see header for description)                               *
 *                                                                                *
 * Authors (alphabetical):                                                        *
 *      Andreas Hoecker <Andreas.Hocker@cern.ch> - CERN, Switzerland              *
 *      Xavier Prudent  <prudent@lapp.in2p3.fr>  - LAPP, France                   *
 *      Helge Voss      <Helge.Voss@cern.ch>     - MPI-K Heidelberg, Germany      *
 *      Kai Voss        <Kai.Voss@cern.ch>       - U. of Victoria, Canada         *
 *                                                                                *
 * Copyright (c) 2005:                                                            *
 *      CERN, Switzerland                                                         * 
 *      U. of Victoria, Canada                                                    * 
 *      MPI-K Heidelberg, Germany                                                 * 
 *      LAPP, Annecy, France                                                      *
 *                                                                                *
 * Redistribution and use in source and binary forms, with or without             *
 * modification, are permitted according to the terms listed in LICENSE           *
 * (http://tmva.sourceforge.net/LICENSE)                                          *
 **********************************************************************************/

//_______________________________________________________________________
//                                                                      
// Begin_Html
/*
  Interface to Clermond-Ferrand artificial neural network 

  <p>
  The CFMlpANN belong to the class of Multilayer Perceptrons (MLP), which are 
  feed-forward networks according to the following propagation schema:<br>
  <center>
  <img vspace=10 src="gif/tmva_mlp.gif" align="bottom" alt="Schema for artificial neural network"> 
  </center>
  The input layer contains as many neurons as input variables used in the MVA.
  The output layer contains two neurons for the signal and background
  event classes. In between the input and output layers are a variable number
  of <i>k</i> hidden layers with arbitrary numbers of neurons. (While the 
  structure of the input and output layers is determined by the problem, the 
  hidden layers can be configured by the user through the option string
  of the method booking.) <br>

  As indicated in the sketch, all neuron inputs to a layer are linear 
  combinations of the neuron output of the previous layer. The transfer
  from input to output within a neuron is performed by means of an "activation
  function". In general, the activation function of a neuron can be 
  zero (deactivated), one (linear), or non-linear. The above example uses
  a sigmoid activation function. The transfer function of the output layer
  is usually linear. As a consequence: an ANN without hidden layer should 
  give identical discrimination power as a linear discriminant analysis (Fisher).
  In case of one hidden layer, the ANN computes a linear combination of 
  sigmoid.  <br>

  The learning method used by the CFMlpANN is only stochastic.
*/
// End_Html
//_______________________________________________________________________

#include <string>
#include <cstdlib>
#include <iostream>

#include "TMatrix.h"
#include "TObjString.h"
#include "Riostream.h"
#include "TMath.h"

#include "TMVA/ClassifierFactory.h"
#include "TMVA/MethodCFMlpANN.h"
#include "TMVA/MethodCFMlpANN_def.h"
#include "TMVA/Tools.h"

REGISTER_METHOD(CFMlpANN)

ClassImp(TMVA::MethodCFMlpANN)

// initialization of global variable
namespace TMVA {
   Int_t MethodCFMlpANN_nsel = 0;
}

TMVA::MethodCFMlpANN* TMVA::MethodCFMlpANN::fgThis = 0;

//_______________________________________________________________________
TMVA::MethodCFMlpANN::MethodCFMlpANN( const TString& jobName,
                                      const TString& methodTitle,
                                      DataSetInfo& theData, 
                                      const TString& theOption,
                                      TDirectory* theTargetDir  ) :
   TMVA::MethodBase( jobName, Types::kCFMlpANN, methodTitle, theData, theOption, theTargetDir  ),
   fNodes(0),
   fYNN(0)
{
   // standard constructor
   // option string: "n_training_cycles:n_hidden_layers"  
   // default is:  n_training_cycles = 5000, n_layers = 4 
   //
   // * note that the number of hidden layers in the NN is: 
   //   n_hidden_layers = n_layers - 2
   //
   // * since there is one input and one output layer. The number of         
   //   nodes (neurons) is predefined to be:
   //   n_nodes[i] = nvars + 1 - i (where i=1..n_layers)                  
   //
   //   with nvars being the number of variables used in the NN.             
   //
   // Hence, the default case is: n_neurons(layer 1 (input)) : nvars       
   //                             n_neurons(layer 2 (hidden)): nvars-1     
   //                             n_neurons(layer 3 (hidden)): nvars-1     
   //                             n_neurons(layer 4 (out))   : 2           
   //
   // This artificial neural network usually needs a relatively large      
   // number of cycles to converge (8000 and more). Overtraining can       
   // be efficienctly tested by comparing the signal and background        
   // output of the NN for the events that were used for training and      
   // an independent data sample (with equal properties). If the separation
   // performance is significantly better for the training sample, the     
   // NN interprets statistical effects, and is hence overtrained. In       
   // this case, the number of cycles should be reduced, or the size       
   // of the training sample increased.                                    
}

//_______________________________________________________________________
TMVA::MethodCFMlpANN::MethodCFMlpANN( DataSetInfo& theData, 
                                      const TString& theWeightFile,  
                                      TDirectory* theTargetDir )
   : TMVA::MethodBase( Types::kCFMlpANN, theData, theWeightFile, theTargetDir ) 
   , fNodes(0)
{
   // constructor from weight file
}

//_______________________________________________________________________
Bool_t TMVA::MethodCFMlpANN::HasAnalysisType( Types::EAnalysisType type, UInt_t numberClasses, UInt_t /*numberTargets*/ )
{
   // CFMlpANN can handle classification with 2 classes 
   if (type == Types::kClassification && numberClasses == 2) return kTRUE;
   return kFALSE;
}

//_______________________________________________________________________
void TMVA::MethodCFMlpANN::DeclareOptions() 
{
   // define the options (their key words) that can be set in the option string 
   // know options: NCycles=xx              :the number of training cycles
   //               HiddenLayser="N-1,N-2"  :the specification of the hidden layers
 
   DeclareOptionRef( fNcycles  =3000,      "NCycles",      "Number of training cycles" );
   DeclareOptionRef( fLayerSpec="N,N-1",   "HiddenLayers", "Specification of hidden layer architecture" );
}

//_______________________________________________________________________
void TMVA::MethodCFMlpANN::ProcessOptions() 
{
   // decode the options in the option string
   fNodes = new Int_t[20]; // number of nodes per layer (maximum 20 layers)
   fNlayers = 2;
   Int_t currentHiddenLayer = 1;
   TString layerSpec(fLayerSpec);
   while(layerSpec.Length()>0) {
      TString sToAdd = "";
      if (layerSpec.First(',')<0) {
         sToAdd = layerSpec;
         layerSpec = "";
      } 
      else {
         sToAdd = layerSpec(0,layerSpec.First(','));
         layerSpec = layerSpec(layerSpec.First(',')+1,layerSpec.Length());
      }
      Int_t nNodes = 0;
      if (sToAdd.BeginsWith("N") || sToAdd.BeginsWith("n")) { sToAdd.Remove(0,1); nNodes = GetNvar(); }
      nNodes += atoi(sToAdd);
      fNodes[currentHiddenLayer++] = nNodes;
      fNlayers++;
   }
   fNodes[0]          = GetNvar(); // number of input nodes
   fNodes[fNlayers-1] = 2;         // number of output nodes

   if (IgnoreEventsWithNegWeightsInTraining()) {
      Log() << kFATAL << "Mechanism to ignore events with negative weights in training not yet available for method: "
            << GetMethodTypeName() 
            << " --> please remove \"IgnoreNegWeightsInTraining\" option from booking string."
            << Endl;
   }

   Log() << kINFO << "Use configuration (nodes per layer): in=";
   for (Int_t i=0; i<fNlayers-1; i++) Log() << kINFO << fNodes[i] << ":";
   Log() << kINFO << fNodes[fNlayers-1] << "=out" << Endl;   

   // some info
   Log() << "Use " << fNcycles << " training cycles" << Endl;

   Int_t nEvtTrain = Data()->GetNTrainingEvents();

   // note that one variable is type
   if (nEvtTrain>0) {
    
      // Data LUT
      fData  = new TMatrix( nEvtTrain, GetNvar() );
      fClass = new vector<Int_t>( nEvtTrain );

      // ---- fill LUTs

      UInt_t ivar;
      for (Int_t ievt=0; ievt<nEvtTrain; ievt++) {
         const Event * ev = GetEvent(ievt);

         // identify signal and background events  
         (*fClass)[ievt] = ev->IsSignal() ? 1 : 2;
      
         // use normalized input Data
         for (ivar=0; ivar<GetNvar(); ivar++) {
            (*fData)( ievt, ivar ) = ev->GetVal(ivar);
         }
      }

      //Log() << kVERBOSE << Data()->GetNEvtSigTrain() << " Signal and " 
      //        << Data()->GetNEvtBkgdTrain() << " background" << " events in trainingTree" << Endl;
   }

}

//_______________________________________________________________________
void TMVA::MethodCFMlpANN::Init( void )
{
   // default initialisation called by all constructors

   // CFMlpANN prefers normalised input variables
   SetNormalised( kTRUE );

   // initialize all pointers
   fgThis = this;  

   // initialize dimensions
   TMVA::MethodCFMlpANN_nsel = 0;  
}

//_______________________________________________________________________
TMVA::MethodCFMlpANN::~MethodCFMlpANN( void )
{
   // destructor
   delete fData;
   delete fClass;
   delete[] fNodes;

   if (fYNN!=0) {
      for (Int_t i=0; i<fNlayers; i++) delete[] fYNN[i];
      delete[] fYNN;
      fYNN=0;
   }
}

//_______________________________________________________________________
void TMVA::MethodCFMlpANN::Train( void )
{
   // training of the Clement-Ferrand NN classifier

   Double_t dumDat(0);
   Int_t ntrain(Data()->GetNTrainingEvents());
   Int_t ntest(0);
   Int_t nvar(GetNvar());
   Int_t nlayers(fNlayers);
   Int_t *nodes = new Int_t[nlayers];
   Int_t ncycles(fNcycles);

   for (Int_t i=0; i<nlayers; i++) nodes[i] = fNodes[i]; // full copy of class member

   if (fYNN != 0) {
      for (Int_t i=0; i<fNlayers; i++) delete[] fYNN[i];
      delete[] fYNN;
      fYNN = 0;
   }
   fYNN = new Double_t*[nlayers];
   for (Int_t layer=0; layer<nlayers; layer++)
      fYNN[layer] = new Double_t[fNodes[layer]];
   
   // please check
#ifndef R__WIN32
   Train_nn( &dumDat, &dumDat, &ntrain, &ntest, &nvar, &nlayers, nodes, &ncycles );
#else
   Log() << kWARNING << "<Train> sorry CFMlpANN does not run on Windows" << Endl;
#endif  

   delete [] nodes;
}

//_______________________________________________________________________
Double_t TMVA::MethodCFMlpANN::GetMvaValue( Double_t* err )
{
   // returns CFMlpANN output (normalised within [0,1])
   Bool_t isOK = kTRUE;

   const Event* ev = GetEvent();

   // copy of input variables
   vector<Double_t> inputVec( GetNvar() );
   for (UInt_t ivar=0; ivar<GetNvar(); ivar++) inputVec[ivar] = ev->GetVal(ivar);

   Double_t myMVA = EvalANN( inputVec, isOK );
   if (!isOK) Log() << kFATAL << "EvalANN returns (!isOK) for event " << Endl;

   // cannot determine error
   if (err != 0) *err = -1;

   return myMVA;
}

//_______________________________________________________________________
Double_t TMVA::MethodCFMlpANN::EvalANN( vector<Double_t>& inVar, Bool_t& isOK )
{
   // evaluates NN value as function of input variables

   // hardcopy of input variables (necessary because they are update later)
   Double_t* xeev = new Double_t[GetNvar()];
   for (UInt_t ivar=0; ivar<GetNvar(); ivar++) xeev[ivar] = inVar[ivar];
  
   // ---- now apply the weights: get NN output
   isOK = kTRUE;
   for (UInt_t jvar=0; jvar<GetNvar(); jvar++) {

      if (fVarn_1.xmax[jvar] < xeev[jvar]) xeev[jvar] = fVarn_1.xmax[jvar];
      if (fVarn_1.xmin[jvar] > xeev[jvar]) xeev[jvar] = fVarn_1.xmin[jvar];
      if (fVarn_1.xmax[jvar] == fVarn_1.xmin[jvar]) {
         isOK = kFALSE;
         xeev[jvar] = 0;
      }
      else {
         xeev[jvar] = xeev[jvar] - ((fVarn_1.xmax[jvar] + fVarn_1.xmin[jvar])/2);    
         xeev[jvar] = xeev[jvar] / ((fVarn_1.xmax[jvar] - fVarn_1.xmin[jvar])/2);    
      }
   }
    
   NN_ava( xeev );

   Double_t retval = 0.5*(1.0 + fYNN[fParam_1.layerm-1][0]);

   delete [] xeev;

   return retval;
}

//_______________________________________________________________________
void  TMVA::MethodCFMlpANN::NN_ava( Double_t* xeev )
{  
   // auxiliary functions
   for (Int_t ivar=0; ivar<fNeur_1.neuron[0]; ivar++) fYNN[0][ivar] = xeev[ivar];

   for (Int_t layer=1; layer<fParam_1.layerm; layer++) {
      for (Int_t j=1; j<=fNeur_1.neuron[layer]; j++) {

         Double_t x = Ww_ref(fNeur_1.ww, layer+1,j); // init with the bias layer

         for (Int_t k=1; k<=fNeur_1.neuron[layer-1]; k++) { // neurons of originating layer
            x += fYNN[layer-1][k-1]*W_ref(fNeur_1.w, layer+1, j, k);
         }
         fYNN[layer][j-1] = NN_fonc( layer, x );
      }
   }  
}

//_______________________________________________________________________
Double_t TMVA::MethodCFMlpANN::NN_fonc( Int_t i, Double_t u ) const
{
   // activation function
   Double_t f(0);
  
   if      (u/fDel_1.temp[i] >  170) f = +1;
   else if (u/fDel_1.temp[i] < -170) f = -1;
   else {
      Double_t yy = TMath::Exp(-u/fDel_1.temp[i]);
      f  = (1 - yy)/(1 + yy);
   }

   return f;
}

//_______________________________________________________________________
void TMVA::MethodCFMlpANN::ReadWeightsFromStream( istream & istr )
{
   // read back the weight from the training from file (stream)
   TString var;

   // read number of variables and classes
   UInt_t nva(0), lclass(0);
   istr >> nva >> lclass;

   if (GetNvar() != nva) // wrong file
      Log() << kFATAL << "<ReadWeightsFromFile> mismatch in number of variables" << Endl;

   // number of output classes must be 2
   if (lclass != 2) // wrong file
      Log() << kFATAL << "<ReadWeightsFromFile> mismatch in number of classes" << Endl;
          
   // check that we are not at the end of the file
   if (istr.eof( ))
      Log() << kFATAL << "<ReadWeightsFromStream> reached EOF prematurely " << Endl;

   // read extrema of input variables
   for (UInt_t ivar=0; ivar<GetNvar(); ivar++) 
      istr >> fVarn_1.xmax[ivar] >> fVarn_1.xmin[ivar];
            
   // read number of layers (sum of: input + output + hidden)
   istr >> fParam_1.layerm;
            
   if (fYNN != 0) {
      for (Int_t i=0; i<fNlayers; i++) delete[] fYNN[i];
      delete[] fYNN;
      fYNN = 0;
   }
   fYNN = new Double_t*[fParam_1.layerm];
   for (Int_t layer=0; layer<fParam_1.layerm; layer++) {              
      // read number of neurons for each layer
      istr >> fNeur_1.neuron[layer];
      fYNN[layer] = new Double_t[fNeur_1.neuron[layer]];
   }
            
   // to read dummy lines
   const Int_t nchar( 100 );
   char* dumchar = new char[nchar];
            
   // read weights
   for (Int_t layer=1; layer<=fParam_1.layerm-1; layer++) {
              
      Int_t nq = fNeur_1.neuron[layer]/10;
      Int_t nr = fNeur_1.neuron[layer] - nq*10;
              
      Int_t kk(0);
      if (nr==0) kk = nq;
      else       kk = nq+1;
              
      for (Int_t k=1; k<=kk; k++) {
         Int_t jmin = 10*k - 9;
         Int_t jmax = 10*k;
         if (fNeur_1.neuron[layer]<jmax) jmax = fNeur_1.neuron[layer];
         for (Int_t j=jmin; j<=jmax; j++) {
            istr >> Ww_ref(fNeur_1.ww, layer+1, j);
         }
         for (Int_t i=1; i<=fNeur_1.neuron[layer-1]; i++) {
            for (Int_t j=jmin; j<=jmax; j++) {
               istr >> W_ref(fNeur_1.w, layer+1, j, i);
            }
         }
         // skip two empty lines
         istr.getline( dumchar, nchar );
      }
   }

   for (Int_t layer=0; layer<fParam_1.layerm; layer++) {
              
      // skip 2 empty lines
      istr.getline( dumchar, nchar );
      istr.getline( dumchar, nchar );
              
      istr >> fDel_1.temp[layer];
   }            

   // sanity check
   if ((Int_t)GetNvar() != fNeur_1.neuron[0]) {
      Log() << kFATAL << "<ReadWeightsFromFile> mismatch in zeroth layer:"
              << GetNvar() << " " << fNeur_1.neuron[0] << Endl;
   }

   fNlayers = fParam_1.layerm;
}

//_______________________________________________________________________
Int_t TMVA::MethodCFMlpANN::DataInterface( Double_t* /*tout2*/, Double_t*  /*tin2*/, 
                                           Int_t* /* icode*/, Int_t*  /*flag*/, 
                                           Int_t*  /*nalire*/, Int_t* nvar, 
                                           Double_t* xpg, Int_t* iclass, Int_t* ikend )
{
   // data interface function 
   
   // icode and ikend are dummies needed to match f2c mlpl3 functions
   *ikend = 0; 

   // retrieve pointer to current object (CFMlpANN must be a singleton class!)
   TMVA::MethodCFMlpANN* opt = TMVA::MethodCFMlpANN::This();

   // sanity checks
   if (0 == xpg) {
      Log() << kFATAL << "ERROR in MethodCFMlpANN_DataInterface zero pointer xpg" << Endl;
   }
   if (*nvar != (Int_t)opt->GetNvar()) {
      Log() << kFATAL << "ERROR in MethodCFMlpANN_DataInterface mismatch in num of variables: "
              << *nvar << " " << opt->GetNvar() << Endl;
   }

   // fill variables
   *iclass = (int)opt->GetClass( TMVA::MethodCFMlpANN_nsel );
   for (UInt_t ivar=0; ivar<opt->GetNvar(); ivar++) 
      xpg[ivar] = (double)opt->GetData( TMVA::MethodCFMlpANN_nsel, ivar );

   ++TMVA::MethodCFMlpANN_nsel;

   return 0;
}

//_______________________________________________________________________
void TMVA::MethodCFMlpANN::WriteWeightsToStream( std::ostream& o ) const
{
   // write number of variables and classes
   o << fParam_1.nvar << "    " << fParam_1.lclass << endl;
   
   // number of output classes must be 2
   if (fParam_1.lclass != 2) // wrong file
      Log() << kFATAL << "<WriteWeightsToStream> mismatch in number of classes" << Endl;


   // write extrema of input variables
   for (Int_t ivar=0; ivar<fParam_1.nvar; ivar++) 
      o << fVarn_1.xmax[ivar] << "   " << fVarn_1.xmin[ivar] << endl;
        
   // write number of layers (sum of: input + output + hidden)
   o << fParam_1.layerm << endl;
        
   for (Int_t layer=0; layer<fParam_1.layerm; layer++) {              
      // write number of neurons for each layer
      o << fNeur_1.neuron[layer] << "     ";
   }
   o << endl;
        
   // write weights
   for (Int_t layer=1; layer<=fParam_1.layerm-1; layer++) { 
          
      Int_t nq = fNeur_1.neuron[layer]/10;
      Int_t nr = fNeur_1.neuron[layer] - nq*10;
          
      Int_t kk(0);
      if (nr==0) kk = nq;
      else       kk = nq+1;
          
      for (Int_t k=1; k<=kk; k++) {
         Int_t jmin = 10*k - 9;
         Int_t jmax = 10*k;
         Int_t i, j;
         if (fNeur_1.neuron[layer]<jmax) jmax = fNeur_1.neuron[layer];
         for (j=jmin; j<=jmax; j++) {
            o << setw(20) << Ww_ref(fNeur_1.ww, layer+1, j) << "  ";
         }
         o << endl;
         for (i=1; i<=fNeur_1.neuron[layer-1]; i++) {
            for (j=jmin; j<=jmax; j++) {
               o << setw(20) << W_ref(fNeur_1.w, layer+1, j, i) << "  ";
            }
            o << endl;
         }
            
         // skip two empty lines
         o << endl << endl;
      }
   }
   for (Int_t layer=0; layer<fParam_1.layerm; layer++) {
          
      // skip 2 empty lines
      o << endl << endl;          
      o << fDel_1.temp[layer] << endl;
   }      
}

//_______________________________________________________________________
void TMVA::MethodCFMlpANN::AddWeightsXMLTo( void* parent ) const 
{
   // write weights to xml file

   void *wght = gTools().xmlengine().NewChild(parent, 0, "Weights");
   gTools().AddAttr(wght,"NVars",fParam_1.nvar);
   gTools().AddAttr(wght,"NClasses",fParam_1.lclass);
   gTools().AddAttr(wght,"NLayers",fParam_1.layerm);
   void* minmaxnode = gTools().xmlengine().NewChild(wght, 0, "VarMinMax");
   stringstream s;
   s.precision( 16 );
   for (Int_t ivar=0; ivar<fParam_1.nvar; ivar++) 
      s << std::scientific << fVarn_1.xmin[ivar] <<  " " << fVarn_1.xmax[ivar] <<  " ";
   gTools().xmlengine().AddRawLine( minmaxnode, s.str().c_str() );
   void* neurons = gTools().xmlengine().NewChild(wght, 0, "NNeurons");
   stringstream n;
   n.precision( 16 );
   for (Int_t layer=0; layer<fParam_1.layerm; layer++)
      n << std::scientific << fNeur_1.neuron[layer] << " ";
   gTools().xmlengine().AddRawLine( neurons, n.str().c_str() );
   for (Int_t layer=1; layer<fParam_1.layerm; layer++) {
      void* layernode = gTools().xmlengine().NewChild(wght, 0, "Layer"+gTools().StringFromInt(layer));
      gTools().AddAttr(layernode,"NNeurons",fNeur_1.neuron[layer]);
      void* neuronnode=NULL;
      for (Int_t neuron=0; neuron<fNeur_1.neuron[layer]; neuron++) {
         neuronnode = gTools().xmlengine().NewChild(layernode,0,"Neuron"+gTools().StringFromInt(neuron));
         stringstream weights;
         weights.precision( 16 );         
         weights << std::scientific << Ww_ref(fNeur_1.ww, layer+1, neuron+1);
         for (Int_t i=0; i<fNeur_1.neuron[layer-1]; i++) {
            weights << " " << std::scientific << W_ref(fNeur_1.w, layer+1, neuron+1, i+1);
         }
         gTools().xmlengine().AddRawLine( neuronnode, weights.str().c_str() );
      }
   }
   void* tempnode = gTools().xmlengine().NewChild(wght, 0, "LayerTemp");
   stringstream temp;
   temp.precision( 16 );
   for (Int_t layer=0; layer<fParam_1.layerm; layer++) {         
       temp << std::scientific << fDel_1.temp[layer] << " ";
   }   
   gTools().xmlengine().AddRawLine(tempnode, temp.str().c_str() );
}
//_______________________________________________________________________
void TMVA::MethodCFMlpANN::ReadWeightsFromXML( void* wghtnode )
{
   // read weights from xml file
   gTools().ReadAttr( wghtnode, "NLayers",fParam_1.layerm );
   void* minmaxnode = gTools().xmlengine().GetChild(wghtnode);
   const char* minmaxcontent = gTools().xmlengine().GetNodeContent(minmaxnode);
   std::stringstream content(minmaxcontent);
   for (UInt_t ivar=0; ivar<GetNvar(); ivar++) 
      content >> fVarn_1.xmin[ivar] >> fVarn_1.xmax[ivar];
   if (fYNN != 0) {
      for (Int_t i=0; i<fNlayers; i++) delete[] fYNN[i];
      delete[] fYNN;
      fYNN = 0;
   }
   fYNN = new Double_t*[fParam_1.layerm];
   void *layernode=gTools().xmlengine().GetNext(minmaxnode);
   const char* neuronscontent = gTools().xmlengine().GetNodeContent(layernode);
   stringstream ncontent(neuronscontent);
   for (Int_t layer=0; layer<fParam_1.layerm; layer++) {              
      // read number of neurons for each layer;
      ncontent >> fNeur_1.neuron[layer];
      fYNN[layer] = new Double_t[fNeur_1.neuron[layer]];
   }
   for (Int_t layer=1; layer<fParam_1.layerm; layer++) {
      layernode=gTools().xmlengine().GetNext(layernode);
      void* neuronnode=NULL;
      neuronnode = gTools().xmlengine().GetChild(layernode);
      for (Int_t neuron=0; neuron<fNeur_1.neuron[layer]; neuron++) {
         const char* neuronweights = gTools().xmlengine().GetNodeContent(neuronnode);
         stringstream weights(neuronweights);
         weights >> Ww_ref(fNeur_1.ww, layer+1, neuron+1);
         for (Int_t i=0; i<fNeur_1.neuron[layer-1]; i++) {
            weights >> W_ref(fNeur_1.w, layer+1, neuron+1, i+1);
         }
         neuronnode=gTools().xmlengine().GetNext(neuronnode);
      }
   } 
   void* tempnode=gTools().xmlengine().GetNext(layernode);
   const char* temp = gTools().xmlengine().GetNodeContent(tempnode);
   stringstream t(temp);
   for (Int_t layer=0; layer<fParam_1.layerm; layer++) {
      t >> fDel_1.temp[layer];
   }
   fNlayers = fParam_1.layerm;
}

//_______________________________________________________________________
void TMVA::MethodCFMlpANN::PrintWeights( std::ostream & o ) const
{
   // write the weights of the neural net

   // write number of variables and classes
   o << "Number of vars " << fParam_1.nvar << endl;
   o << "Output nodes   " << fParam_1.lclass << endl;
   
   // write extrema of input variables
   for (Int_t ivar=0; ivar<fParam_1.nvar; ivar++) 
      o << "Var " << ivar << " [" << fVarn_1.xmin[ivar] << " - " << fVarn_1.xmax[ivar] << "]" << endl;
        
   // write number of layers (sum of: input + output + hidden)
   o << "Number of layers " << fParam_1.layerm << endl;
   
   o << "Nodes per layer ";
   for (Int_t layer=0; layer<fParam_1.layerm; layer++)
      // write number of neurons for each layer
      o << fNeur_1.neuron[layer] << "     ";   
   o << endl;
        
   // write weights
   for (Int_t layer=1; layer<=fParam_1.layerm-1; layer++) { 
          
      Int_t nq = fNeur_1.neuron[layer]/10;
      Int_t nr = fNeur_1.neuron[layer] - nq*10;
          
      Int_t kk(0);
      if (nr==0) kk = nq;
      else       kk = nq+1;
          
      for (Int_t k=1; k<=kk; k++) {
         Int_t jmin = 10*k - 9;
         Int_t jmax = 10*k;
         Int_t i, j;
         if (fNeur_1.neuron[layer]<jmax) jmax = fNeur_1.neuron[layer];
         for (j=jmin; j<=jmax; j++) {

            //o << fNeur_1.ww[j*max_nLayers_ + layer - 6] << "   ";
            o << Ww_ref(fNeur_1.ww, layer+1, j) << "   ";

         }
         o << endl;
         //for (i=1; i<=fNeur_1.neuron[layer-1]; i++) {
         for (i=1; i<=fNeur_1.neuron[layer-1]; i++) {
            for (j=jmin; j<=jmax; j++) {
               //               o << fNeur_1.w[(i*max_nNodes_ + j)*max_nLayers_ + layer - 186] << "   ";
               o << W_ref(fNeur_1.w, layer+1, j, i) << "   ";
            }
            o << endl;
         }
            
         // skip two empty lines
         o << endl;
      }
   }
   for (Int_t layer=0; layer<fParam_1.layerm; layer++) {
      o << "Del.temp in layer " << layer << " :  " << fDel_1.temp[layer] << endl;
   }      
}

//_______________________________________________________________________
void TMVA::MethodCFMlpANN::MakeClassSpecific( std::ostream& fout, const TString& className ) const
{
   // write specific classifier response
   fout << "   // not implemented for class: \"" << className << "\"" << endl;
   fout << "};" << endl;
}

//_______________________________________________________________________
void TMVA::MethodCFMlpANN::MakeClassSpecificHeader( std::ostream& , const TString&  ) const
{
   // write specific classifier response for header
}

//_______________________________________________________________________
void TMVA::MethodCFMlpANN::GetHelpMessage() const
{
   // get help message text
   //
   // typical length of text line: 
   //         "|--------------------------------------------------------------|"
   Log() << Endl;
   Log() << gTools().Color("bold") << "--- Short description:" << gTools().Color("reset") << Endl;
   Log() << Endl;
   Log() << "<None>" << Endl;
   Log() << Endl;
   Log() << gTools().Color("bold") << "--- Performance optimisation:" << gTools().Color("reset") << Endl;
   Log() << Endl;
   Log() << "<None>" << Endl;
   Log() << Endl;
   Log() << gTools().Color("bold") << "--- Performance tuning via configuration options:" << gTools().Color("reset") << Endl;
   Log() << Endl;
   Log() << "<None>" << Endl;
}
 MethodCFMlpANN.cxx:1
 MethodCFMlpANN.cxx:2
 MethodCFMlpANN.cxx:3
 MethodCFMlpANN.cxx:4
 MethodCFMlpANN.cxx:5
 MethodCFMlpANN.cxx:6
 MethodCFMlpANN.cxx:7
 MethodCFMlpANN.cxx:8
 MethodCFMlpANN.cxx:9
 MethodCFMlpANN.cxx:10
 MethodCFMlpANN.cxx:11
 MethodCFMlpANN.cxx:12
 MethodCFMlpANN.cxx:13
 MethodCFMlpANN.cxx:14
 MethodCFMlpANN.cxx:15
 MethodCFMlpANN.cxx:16
 MethodCFMlpANN.cxx:17
 MethodCFMlpANN.cxx:18
 MethodCFMlpANN.cxx:19
 MethodCFMlpANN.cxx:20
 MethodCFMlpANN.cxx:21
 MethodCFMlpANN.cxx:22
 MethodCFMlpANN.cxx:23
 MethodCFMlpANN.cxx:24
 MethodCFMlpANN.cxx:25
 MethodCFMlpANN.cxx:26
 MethodCFMlpANN.cxx:27
 MethodCFMlpANN.cxx:28
 MethodCFMlpANN.cxx:29
 MethodCFMlpANN.cxx:30
 MethodCFMlpANN.cxx:31
 MethodCFMlpANN.cxx:32
 MethodCFMlpANN.cxx:33
 MethodCFMlpANN.cxx:34
 MethodCFMlpANN.cxx:35
 MethodCFMlpANN.cxx:36
 MethodCFMlpANN.cxx:37
 MethodCFMlpANN.cxx:38
 MethodCFMlpANN.cxx:39
 MethodCFMlpANN.cxx:40
 MethodCFMlpANN.cxx:41
 MethodCFMlpANN.cxx:42
 MethodCFMlpANN.cxx:43
 MethodCFMlpANN.cxx:44
 MethodCFMlpANN.cxx:45
 MethodCFMlpANN.cxx:46
 MethodCFMlpANN.cxx:47
 MethodCFMlpANN.cxx:48
 MethodCFMlpANN.cxx:49
 MethodCFMlpANN.cxx:50
 MethodCFMlpANN.cxx:51
 MethodCFMlpANN.cxx:52
 MethodCFMlpANN.cxx:53
 MethodCFMlpANN.cxx:54
 MethodCFMlpANN.cxx:55
 MethodCFMlpANN.cxx:56
 MethodCFMlpANN.cxx:57
 MethodCFMlpANN.cxx:58
 MethodCFMlpANN.cxx:59
 MethodCFMlpANN.cxx:60
 MethodCFMlpANN.cxx:61
 MethodCFMlpANN.cxx:62
 MethodCFMlpANN.cxx:63
 MethodCFMlpANN.cxx:64
 MethodCFMlpANN.cxx:65
 MethodCFMlpANN.cxx:66
 MethodCFMlpANN.cxx:67
 MethodCFMlpANN.cxx:68
 MethodCFMlpANN.cxx:69
 MethodCFMlpANN.cxx:70
 MethodCFMlpANN.cxx:71
 MethodCFMlpANN.cxx:72
 MethodCFMlpANN.cxx:73
 MethodCFMlpANN.cxx:74
 MethodCFMlpANN.cxx:75
 MethodCFMlpANN.cxx:76
 MethodCFMlpANN.cxx:77
 MethodCFMlpANN.cxx:78
 MethodCFMlpANN.cxx:79
 MethodCFMlpANN.cxx:80
 MethodCFMlpANN.cxx:81
 MethodCFMlpANN.cxx:82
 MethodCFMlpANN.cxx:83
 MethodCFMlpANN.cxx:84
 MethodCFMlpANN.cxx:85
 MethodCFMlpANN.cxx:86
 MethodCFMlpANN.cxx:87
 MethodCFMlpANN.cxx:88
 MethodCFMlpANN.cxx:89
 MethodCFMlpANN.cxx:90
 MethodCFMlpANN.cxx:91
 MethodCFMlpANN.cxx:92
 MethodCFMlpANN.cxx:93
 MethodCFMlpANN.cxx:94
 MethodCFMlpANN.cxx:95
 MethodCFMlpANN.cxx:96
 MethodCFMlpANN.cxx:97
 MethodCFMlpANN.cxx:98
 MethodCFMlpANN.cxx:99
 MethodCFMlpANN.cxx:100
 MethodCFMlpANN.cxx:101
 MethodCFMlpANN.cxx:102
 MethodCFMlpANN.cxx:103
 MethodCFMlpANN.cxx:104
 MethodCFMlpANN.cxx:105
 MethodCFMlpANN.cxx:106
 MethodCFMlpANN.cxx:107
 MethodCFMlpANN.cxx:108
 MethodCFMlpANN.cxx:109
 MethodCFMlpANN.cxx:110
 MethodCFMlpANN.cxx:111
 MethodCFMlpANN.cxx:112
 MethodCFMlpANN.cxx:113
 MethodCFMlpANN.cxx:114
 MethodCFMlpANN.cxx:115
 MethodCFMlpANN.cxx:116
 MethodCFMlpANN.cxx:117
 MethodCFMlpANN.cxx:118
 MethodCFMlpANN.cxx:119
 MethodCFMlpANN.cxx:120
 MethodCFMlpANN.cxx:121
 MethodCFMlpANN.cxx:122
 MethodCFMlpANN.cxx:123
 MethodCFMlpANN.cxx:124
 MethodCFMlpANN.cxx:125
 MethodCFMlpANN.cxx:126
 MethodCFMlpANN.cxx:127
 MethodCFMlpANN.cxx:128
 MethodCFMlpANN.cxx:129
 MethodCFMlpANN.cxx:130
 MethodCFMlpANN.cxx:131
 MethodCFMlpANN.cxx:132
 MethodCFMlpANN.cxx:133
 MethodCFMlpANN.cxx:134
 MethodCFMlpANN.cxx:135
 MethodCFMlpANN.cxx:136
 MethodCFMlpANN.cxx:137
 MethodCFMlpANN.cxx:138
 MethodCFMlpANN.cxx:139
 MethodCFMlpANN.cxx:140
 MethodCFMlpANN.cxx:141
 MethodCFMlpANN.cxx:142
 MethodCFMlpANN.cxx:143
 MethodCFMlpANN.cxx:144
 MethodCFMlpANN.cxx:145
 MethodCFMlpANN.cxx:146
 MethodCFMlpANN.cxx:147
 MethodCFMlpANN.cxx:148
 MethodCFMlpANN.cxx:149
 MethodCFMlpANN.cxx:150
 MethodCFMlpANN.cxx:151
 MethodCFMlpANN.cxx:152
 MethodCFMlpANN.cxx:153
 MethodCFMlpANN.cxx:154
 MethodCFMlpANN.cxx:155
 MethodCFMlpANN.cxx:156
 MethodCFMlpANN.cxx:157
 MethodCFMlpANN.cxx:158
 MethodCFMlpANN.cxx:159
 MethodCFMlpANN.cxx:160
 MethodCFMlpANN.cxx:161
 MethodCFMlpANN.cxx:162
 MethodCFMlpANN.cxx:163
 MethodCFMlpANN.cxx:164
 MethodCFMlpANN.cxx:165
 MethodCFMlpANN.cxx:166
 MethodCFMlpANN.cxx:167
 MethodCFMlpANN.cxx:168
 MethodCFMlpANN.cxx:169
 MethodCFMlpANN.cxx:170
 MethodCFMlpANN.cxx:171
 MethodCFMlpANN.cxx:172
 MethodCFMlpANN.cxx:173
 MethodCFMlpANN.cxx:174
 MethodCFMlpANN.cxx:175
 MethodCFMlpANN.cxx:176
 MethodCFMlpANN.cxx:177
 MethodCFMlpANN.cxx:178
 MethodCFMlpANN.cxx:179
 MethodCFMlpANN.cxx:180
 MethodCFMlpANN.cxx:181
 MethodCFMlpANN.cxx:182
 MethodCFMlpANN.cxx:183
 MethodCFMlpANN.cxx:184
 MethodCFMlpANN.cxx:185
 MethodCFMlpANN.cxx:186
 MethodCFMlpANN.cxx:187
 MethodCFMlpANN.cxx:188
 MethodCFMlpANN.cxx:189
 MethodCFMlpANN.cxx:190
 MethodCFMlpANN.cxx:191
 MethodCFMlpANN.cxx:192
 MethodCFMlpANN.cxx:193
 MethodCFMlpANN.cxx:194
 MethodCFMlpANN.cxx:195
 MethodCFMlpANN.cxx:196
 MethodCFMlpANN.cxx:197
 MethodCFMlpANN.cxx:198
 MethodCFMlpANN.cxx:199
 MethodCFMlpANN.cxx:200
 MethodCFMlpANN.cxx:201
 MethodCFMlpANN.cxx:202
 MethodCFMlpANN.cxx:203
 MethodCFMlpANN.cxx:204
 MethodCFMlpANN.cxx:205
 MethodCFMlpANN.cxx:206
 MethodCFMlpANN.cxx:207
 MethodCFMlpANN.cxx:208
 MethodCFMlpANN.cxx:209
 MethodCFMlpANN.cxx:210
 MethodCFMlpANN.cxx:211
 MethodCFMlpANN.cxx:212
 MethodCFMlpANN.cxx:213
 MethodCFMlpANN.cxx:214
 MethodCFMlpANN.cxx:215
 MethodCFMlpANN.cxx:216
 MethodCFMlpANN.cxx:217
 MethodCFMlpANN.cxx:218
 MethodCFMlpANN.cxx:219
 MethodCFMlpANN.cxx:220
 MethodCFMlpANN.cxx:221
 MethodCFMlpANN.cxx:222
 MethodCFMlpANN.cxx:223
 MethodCFMlpANN.cxx:224
 MethodCFMlpANN.cxx:225
 MethodCFMlpANN.cxx:226
 MethodCFMlpANN.cxx:227
 MethodCFMlpANN.cxx:228
 MethodCFMlpANN.cxx:229
 MethodCFMlpANN.cxx:230
 MethodCFMlpANN.cxx:231
 MethodCFMlpANN.cxx:232
 MethodCFMlpANN.cxx:233
 MethodCFMlpANN.cxx:234
 MethodCFMlpANN.cxx:235
 MethodCFMlpANN.cxx:236
 MethodCFMlpANN.cxx:237
 MethodCFMlpANN.cxx:238
 MethodCFMlpANN.cxx:239
 MethodCFMlpANN.cxx:240
 MethodCFMlpANN.cxx:241
 MethodCFMlpANN.cxx:242
 MethodCFMlpANN.cxx:243
 MethodCFMlpANN.cxx:244
 MethodCFMlpANN.cxx:245
 MethodCFMlpANN.cxx:246
 MethodCFMlpANN.cxx:247
 MethodCFMlpANN.cxx:248
 MethodCFMlpANN.cxx:249
 MethodCFMlpANN.cxx:250
 MethodCFMlpANN.cxx:251
 MethodCFMlpANN.cxx:252
 MethodCFMlpANN.cxx:253
 MethodCFMlpANN.cxx:254
 MethodCFMlpANN.cxx:255
 MethodCFMlpANN.cxx:256
 MethodCFMlpANN.cxx:257
 MethodCFMlpANN.cxx:258
 MethodCFMlpANN.cxx:259
 MethodCFMlpANN.cxx:260
 MethodCFMlpANN.cxx:261
 MethodCFMlpANN.cxx:262
 MethodCFMlpANN.cxx:263
 MethodCFMlpANN.cxx:264
 MethodCFMlpANN.cxx:265
 MethodCFMlpANN.cxx:266
 MethodCFMlpANN.cxx:267
 MethodCFMlpANN.cxx:268
 MethodCFMlpANN.cxx:269
 MethodCFMlpANN.cxx:270
 MethodCFMlpANN.cxx:271
 MethodCFMlpANN.cxx:272
 MethodCFMlpANN.cxx:273
 MethodCFMlpANN.cxx:274
 MethodCFMlpANN.cxx:275
 MethodCFMlpANN.cxx:276
 MethodCFMlpANN.cxx:277
 MethodCFMlpANN.cxx:278
 MethodCFMlpANN.cxx:279
 MethodCFMlpANN.cxx:280
 MethodCFMlpANN.cxx:281
 MethodCFMlpANN.cxx:282
 MethodCFMlpANN.cxx:283
 MethodCFMlpANN.cxx:284
 MethodCFMlpANN.cxx:285
 MethodCFMlpANN.cxx:286
 MethodCFMlpANN.cxx:287
 MethodCFMlpANN.cxx:288
 MethodCFMlpANN.cxx:289
 MethodCFMlpANN.cxx:290
 MethodCFMlpANN.cxx:291
 MethodCFMlpANN.cxx:292
 MethodCFMlpANN.cxx:293
 MethodCFMlpANN.cxx:294
 MethodCFMlpANN.cxx:295
 MethodCFMlpANN.cxx:296
 MethodCFMlpANN.cxx:297
 MethodCFMlpANN.cxx:298
 MethodCFMlpANN.cxx:299
 MethodCFMlpANN.cxx:300
 MethodCFMlpANN.cxx:301
 MethodCFMlpANN.cxx:302
 MethodCFMlpANN.cxx:303
 MethodCFMlpANN.cxx:304
 MethodCFMlpANN.cxx:305
 MethodCFMlpANN.cxx:306
 MethodCFMlpANN.cxx:307
 MethodCFMlpANN.cxx:308
 MethodCFMlpANN.cxx:309
 MethodCFMlpANN.cxx:310
 MethodCFMlpANN.cxx:311
 MethodCFMlpANN.cxx:312
 MethodCFMlpANN.cxx:313
 MethodCFMlpANN.cxx:314
 MethodCFMlpANN.cxx:315
 MethodCFMlpANN.cxx:316
 MethodCFMlpANN.cxx:317
 MethodCFMlpANN.cxx:318
 MethodCFMlpANN.cxx:319
 MethodCFMlpANN.cxx:320
 MethodCFMlpANN.cxx:321
 MethodCFMlpANN.cxx:322
 MethodCFMlpANN.cxx:323
 MethodCFMlpANN.cxx:324
 MethodCFMlpANN.cxx:325
 MethodCFMlpANN.cxx:326
 MethodCFMlpANN.cxx:327
 MethodCFMlpANN.cxx:328
 MethodCFMlpANN.cxx:329
 MethodCFMlpANN.cxx:330
 MethodCFMlpANN.cxx:331
 MethodCFMlpANN.cxx:332
 MethodCFMlpANN.cxx:333
 MethodCFMlpANN.cxx:334
 MethodCFMlpANN.cxx:335
 MethodCFMlpANN.cxx:336
 MethodCFMlpANN.cxx:337
 MethodCFMlpANN.cxx:338
 MethodCFMlpANN.cxx:339
 MethodCFMlpANN.cxx:340
 MethodCFMlpANN.cxx:341
 MethodCFMlpANN.cxx:342
 MethodCFMlpANN.cxx:343
 MethodCFMlpANN.cxx:344
 MethodCFMlpANN.cxx:345
 MethodCFMlpANN.cxx:346
 MethodCFMlpANN.cxx:347
 MethodCFMlpANN.cxx:348
 MethodCFMlpANN.cxx:349
 MethodCFMlpANN.cxx:350
 MethodCFMlpANN.cxx:351
 MethodCFMlpANN.cxx:352
 MethodCFMlpANN.cxx:353
 MethodCFMlpANN.cxx:354
 MethodCFMlpANN.cxx:355
 MethodCFMlpANN.cxx:356
 MethodCFMlpANN.cxx:357
 MethodCFMlpANN.cxx:358
 MethodCFMlpANN.cxx:359
 MethodCFMlpANN.cxx:360
 MethodCFMlpANN.cxx:361
 MethodCFMlpANN.cxx:362
 MethodCFMlpANN.cxx:363
 MethodCFMlpANN.cxx:364
 MethodCFMlpANN.cxx:365
 MethodCFMlpANN.cxx:366
 MethodCFMlpANN.cxx:367
 MethodCFMlpANN.cxx:368
 MethodCFMlpANN.cxx:369
 MethodCFMlpANN.cxx:370
 MethodCFMlpANN.cxx:371
 MethodCFMlpANN.cxx:372
 MethodCFMlpANN.cxx:373
 MethodCFMlpANN.cxx:374
 MethodCFMlpANN.cxx:375
 MethodCFMlpANN.cxx:376
 MethodCFMlpANN.cxx:377
 MethodCFMlpANN.cxx:378
 MethodCFMlpANN.cxx:379
 MethodCFMlpANN.cxx:380
 MethodCFMlpANN.cxx:381
 MethodCFMlpANN.cxx:382
 MethodCFMlpANN.cxx:383
 MethodCFMlpANN.cxx:384
 MethodCFMlpANN.cxx:385
 MethodCFMlpANN.cxx:386
 MethodCFMlpANN.cxx:387
 MethodCFMlpANN.cxx:388
 MethodCFMlpANN.cxx:389
 MethodCFMlpANN.cxx:390
 MethodCFMlpANN.cxx:391
 MethodCFMlpANN.cxx:392
 MethodCFMlpANN.cxx:393
 MethodCFMlpANN.cxx:394
 MethodCFMlpANN.cxx:395
 MethodCFMlpANN.cxx:396
 MethodCFMlpANN.cxx:397
 MethodCFMlpANN.cxx:398
 MethodCFMlpANN.cxx:399
 MethodCFMlpANN.cxx:400
 MethodCFMlpANN.cxx:401
 MethodCFMlpANN.cxx:402
 MethodCFMlpANN.cxx:403
 MethodCFMlpANN.cxx:404
 MethodCFMlpANN.cxx:405
 MethodCFMlpANN.cxx:406
 MethodCFMlpANN.cxx:407
 MethodCFMlpANN.cxx:408
 MethodCFMlpANN.cxx:409
 MethodCFMlpANN.cxx:410
 MethodCFMlpANN.cxx:411
 MethodCFMlpANN.cxx:412
 MethodCFMlpANN.cxx:413
 MethodCFMlpANN.cxx:414
 MethodCFMlpANN.cxx:415
 MethodCFMlpANN.cxx:416
 MethodCFMlpANN.cxx:417
 MethodCFMlpANN.cxx:418
 MethodCFMlpANN.cxx:419
 MethodCFMlpANN.cxx:420
 MethodCFMlpANN.cxx:421
 MethodCFMlpANN.cxx:422
 MethodCFMlpANN.cxx:423
 MethodCFMlpANN.cxx:424
 MethodCFMlpANN.cxx:425
 MethodCFMlpANN.cxx:426
 MethodCFMlpANN.cxx:427
 MethodCFMlpANN.cxx:428
 MethodCFMlpANN.cxx:429
 MethodCFMlpANN.cxx:430
 MethodCFMlpANN.cxx:431
 MethodCFMlpANN.cxx:432
 MethodCFMlpANN.cxx:433
 MethodCFMlpANN.cxx:434
 MethodCFMlpANN.cxx:435
 MethodCFMlpANN.cxx:436
 MethodCFMlpANN.cxx:437
 MethodCFMlpANN.cxx:438
 MethodCFMlpANN.cxx:439
 MethodCFMlpANN.cxx:440
 MethodCFMlpANN.cxx:441
 MethodCFMlpANN.cxx:442
 MethodCFMlpANN.cxx:443
 MethodCFMlpANN.cxx:444
 MethodCFMlpANN.cxx:445
 MethodCFMlpANN.cxx:446
 MethodCFMlpANN.cxx:447
 MethodCFMlpANN.cxx:448
 MethodCFMlpANN.cxx:449
 MethodCFMlpANN.cxx:450
 MethodCFMlpANN.cxx:451
 MethodCFMlpANN.cxx:452
 MethodCFMlpANN.cxx:453
 MethodCFMlpANN.cxx:454
 MethodCFMlpANN.cxx:455
 MethodCFMlpANN.cxx:456
 MethodCFMlpANN.cxx:457
 MethodCFMlpANN.cxx:458
 MethodCFMlpANN.cxx:459
 MethodCFMlpANN.cxx:460
 MethodCFMlpANN.cxx:461
 MethodCFMlpANN.cxx:462
 MethodCFMlpANN.cxx:463
 MethodCFMlpANN.cxx:464
 MethodCFMlpANN.cxx:465
 MethodCFMlpANN.cxx:466
 MethodCFMlpANN.cxx:467
 MethodCFMlpANN.cxx:468
 MethodCFMlpANN.cxx:469
 MethodCFMlpANN.cxx:470
 MethodCFMlpANN.cxx:471
 MethodCFMlpANN.cxx:472
 MethodCFMlpANN.cxx:473
 MethodCFMlpANN.cxx:474
 MethodCFMlpANN.cxx:475
 MethodCFMlpANN.cxx:476
 MethodCFMlpANN.cxx:477
 MethodCFMlpANN.cxx:478
 MethodCFMlpANN.cxx:479
 MethodCFMlpANN.cxx:480
 MethodCFMlpANN.cxx:481
 MethodCFMlpANN.cxx:482
 MethodCFMlpANN.cxx:483
 MethodCFMlpANN.cxx:484
 MethodCFMlpANN.cxx:485
 MethodCFMlpANN.cxx:486
 MethodCFMlpANN.cxx:487
 MethodCFMlpANN.cxx:488
 MethodCFMlpANN.cxx:489
 MethodCFMlpANN.cxx:490
 MethodCFMlpANN.cxx:491
 MethodCFMlpANN.cxx:492
 MethodCFMlpANN.cxx:493
 MethodCFMlpANN.cxx:494
 MethodCFMlpANN.cxx:495
 MethodCFMlpANN.cxx:496
 MethodCFMlpANN.cxx:497
 MethodCFMlpANN.cxx:498
 MethodCFMlpANN.cxx:499
 MethodCFMlpANN.cxx:500
 MethodCFMlpANN.cxx:501
 MethodCFMlpANN.cxx:502
 MethodCFMlpANN.cxx:503
 MethodCFMlpANN.cxx:504
 MethodCFMlpANN.cxx:505
 MethodCFMlpANN.cxx:506
 MethodCFMlpANN.cxx:507
 MethodCFMlpANN.cxx:508
 MethodCFMlpANN.cxx:509
 MethodCFMlpANN.cxx:510
 MethodCFMlpANN.cxx:511
 MethodCFMlpANN.cxx:512
 MethodCFMlpANN.cxx:513
 MethodCFMlpANN.cxx:514
 MethodCFMlpANN.cxx:515
 MethodCFMlpANN.cxx:516
 MethodCFMlpANN.cxx:517
 MethodCFMlpANN.cxx:518
 MethodCFMlpANN.cxx:519
 MethodCFMlpANN.cxx:520
 MethodCFMlpANN.cxx:521
 MethodCFMlpANN.cxx:522
 MethodCFMlpANN.cxx:523
 MethodCFMlpANN.cxx:524
 MethodCFMlpANN.cxx:525
 MethodCFMlpANN.cxx:526
 MethodCFMlpANN.cxx:527
 MethodCFMlpANN.cxx:528
 MethodCFMlpANN.cxx:529
 MethodCFMlpANN.cxx:530
 MethodCFMlpANN.cxx:531
 MethodCFMlpANN.cxx:532
 MethodCFMlpANN.cxx:533
 MethodCFMlpANN.cxx:534
 MethodCFMlpANN.cxx:535
 MethodCFMlpANN.cxx:536
 MethodCFMlpANN.cxx:537
 MethodCFMlpANN.cxx:538
 MethodCFMlpANN.cxx:539
 MethodCFMlpANN.cxx:540
 MethodCFMlpANN.cxx:541
 MethodCFMlpANN.cxx:542
 MethodCFMlpANN.cxx:543
 MethodCFMlpANN.cxx:544
 MethodCFMlpANN.cxx:545
 MethodCFMlpANN.cxx:546
 MethodCFMlpANN.cxx:547
 MethodCFMlpANN.cxx:548
 MethodCFMlpANN.cxx:549
 MethodCFMlpANN.cxx:550
 MethodCFMlpANN.cxx:551
 MethodCFMlpANN.cxx:552
 MethodCFMlpANN.cxx:553
 MethodCFMlpANN.cxx:554
 MethodCFMlpANN.cxx:555
 MethodCFMlpANN.cxx:556
 MethodCFMlpANN.cxx:557
 MethodCFMlpANN.cxx:558
 MethodCFMlpANN.cxx:559
 MethodCFMlpANN.cxx:560
 MethodCFMlpANN.cxx:561
 MethodCFMlpANN.cxx:562
 MethodCFMlpANN.cxx:563
 MethodCFMlpANN.cxx:564
 MethodCFMlpANN.cxx:565
 MethodCFMlpANN.cxx:566
 MethodCFMlpANN.cxx:567
 MethodCFMlpANN.cxx:568
 MethodCFMlpANN.cxx:569
 MethodCFMlpANN.cxx:570
 MethodCFMlpANN.cxx:571
 MethodCFMlpANN.cxx:572
 MethodCFMlpANN.cxx:573
 MethodCFMlpANN.cxx:574
 MethodCFMlpANN.cxx:575
 MethodCFMlpANN.cxx:576
 MethodCFMlpANN.cxx:577
 MethodCFMlpANN.cxx:578
 MethodCFMlpANN.cxx:579
 MethodCFMlpANN.cxx:580
 MethodCFMlpANN.cxx:581
 MethodCFMlpANN.cxx:582
 MethodCFMlpANN.cxx:583
 MethodCFMlpANN.cxx:584
 MethodCFMlpANN.cxx:585
 MethodCFMlpANN.cxx:586
 MethodCFMlpANN.cxx:587
 MethodCFMlpANN.cxx:588
 MethodCFMlpANN.cxx:589
 MethodCFMlpANN.cxx:590
 MethodCFMlpANN.cxx:591
 MethodCFMlpANN.cxx:592
 MethodCFMlpANN.cxx:593
 MethodCFMlpANN.cxx:594
 MethodCFMlpANN.cxx:595
 MethodCFMlpANN.cxx:596
 MethodCFMlpANN.cxx:597
 MethodCFMlpANN.cxx:598
 MethodCFMlpANN.cxx:599
 MethodCFMlpANN.cxx:600
 MethodCFMlpANN.cxx:601
 MethodCFMlpANN.cxx:602
 MethodCFMlpANN.cxx:603
 MethodCFMlpANN.cxx:604
 MethodCFMlpANN.cxx:605
 MethodCFMlpANN.cxx:606
 MethodCFMlpANN.cxx:607
 MethodCFMlpANN.cxx:608
 MethodCFMlpANN.cxx:609
 MethodCFMlpANN.cxx:610
 MethodCFMlpANN.cxx:611
 MethodCFMlpANN.cxx:612
 MethodCFMlpANN.cxx:613
 MethodCFMlpANN.cxx:614
 MethodCFMlpANN.cxx:615
 MethodCFMlpANN.cxx:616
 MethodCFMlpANN.cxx:617
 MethodCFMlpANN.cxx:618
 MethodCFMlpANN.cxx:619
 MethodCFMlpANN.cxx:620
 MethodCFMlpANN.cxx:621
 MethodCFMlpANN.cxx:622
 MethodCFMlpANN.cxx:623
 MethodCFMlpANN.cxx:624
 MethodCFMlpANN.cxx:625
 MethodCFMlpANN.cxx:626
 MethodCFMlpANN.cxx:627
 MethodCFMlpANN.cxx:628
 MethodCFMlpANN.cxx:629
 MethodCFMlpANN.cxx:630
 MethodCFMlpANN.cxx:631
 MethodCFMlpANN.cxx:632
 MethodCFMlpANN.cxx:633
 MethodCFMlpANN.cxx:634
 MethodCFMlpANN.cxx:635
 MethodCFMlpANN.cxx:636
 MethodCFMlpANN.cxx:637
 MethodCFMlpANN.cxx:638
 MethodCFMlpANN.cxx:639
 MethodCFMlpANN.cxx:640
 MethodCFMlpANN.cxx:641
 MethodCFMlpANN.cxx:642
 MethodCFMlpANN.cxx:643
 MethodCFMlpANN.cxx:644
 MethodCFMlpANN.cxx:645
 MethodCFMlpANN.cxx:646
 MethodCFMlpANN.cxx:647
 MethodCFMlpANN.cxx:648
 MethodCFMlpANN.cxx:649
 MethodCFMlpANN.cxx:650
 MethodCFMlpANN.cxx:651
 MethodCFMlpANN.cxx:652
 MethodCFMlpANN.cxx:653
 MethodCFMlpANN.cxx:654
 MethodCFMlpANN.cxx:655
 MethodCFMlpANN.cxx:656
 MethodCFMlpANN.cxx:657
 MethodCFMlpANN.cxx:658
 MethodCFMlpANN.cxx:659
 MethodCFMlpANN.cxx:660
 MethodCFMlpANN.cxx:661
 MethodCFMlpANN.cxx:662
 MethodCFMlpANN.cxx:663
 MethodCFMlpANN.cxx:664
 MethodCFMlpANN.cxx:665
 MethodCFMlpANN.cxx:666
 MethodCFMlpANN.cxx:667
 MethodCFMlpANN.cxx:668
 MethodCFMlpANN.cxx:669
 MethodCFMlpANN.cxx:670
 MethodCFMlpANN.cxx:671
 MethodCFMlpANN.cxx:672
 MethodCFMlpANN.cxx:673
 MethodCFMlpANN.cxx:674
 MethodCFMlpANN.cxx:675
 MethodCFMlpANN.cxx:676
 MethodCFMlpANN.cxx:677
 MethodCFMlpANN.cxx:678
 MethodCFMlpANN.cxx:679
 MethodCFMlpANN.cxx:680
 MethodCFMlpANN.cxx:681
 MethodCFMlpANN.cxx:682
 MethodCFMlpANN.cxx:683
 MethodCFMlpANN.cxx:684
 MethodCFMlpANN.cxx:685
 MethodCFMlpANN.cxx:686
 MethodCFMlpANN.cxx:687
 MethodCFMlpANN.cxx:688
 MethodCFMlpANN.cxx:689
 MethodCFMlpANN.cxx:690
 MethodCFMlpANN.cxx:691
 MethodCFMlpANN.cxx:692
 MethodCFMlpANN.cxx:693
 MethodCFMlpANN.cxx:694
 MethodCFMlpANN.cxx:695
 MethodCFMlpANN.cxx:696
 MethodCFMlpANN.cxx:697
 MethodCFMlpANN.cxx:698
 MethodCFMlpANN.cxx:699
 MethodCFMlpANN.cxx:700
 MethodCFMlpANN.cxx:701
 MethodCFMlpANN.cxx:702
 MethodCFMlpANN.cxx:703
 MethodCFMlpANN.cxx:704
 MethodCFMlpANN.cxx:705
 MethodCFMlpANN.cxx:706
 MethodCFMlpANN.cxx:707
 MethodCFMlpANN.cxx:708
 MethodCFMlpANN.cxx:709
 MethodCFMlpANN.cxx:710
 MethodCFMlpANN.cxx:711
 MethodCFMlpANN.cxx:712
 MethodCFMlpANN.cxx:713
 MethodCFMlpANN.cxx:714
 MethodCFMlpANN.cxx:715
 MethodCFMlpANN.cxx:716
 MethodCFMlpANN.cxx:717
 MethodCFMlpANN.cxx:718
 MethodCFMlpANN.cxx:719
 MethodCFMlpANN.cxx:720
 MethodCFMlpANN.cxx:721
 MethodCFMlpANN.cxx:722
 MethodCFMlpANN.cxx:723
 MethodCFMlpANN.cxx:724
 MethodCFMlpANN.cxx:725
 MethodCFMlpANN.cxx:726
 MethodCFMlpANN.cxx:727
 MethodCFMlpANN.cxx:728
 MethodCFMlpANN.cxx:729
 MethodCFMlpANN.cxx:730
 MethodCFMlpANN.cxx:731
 MethodCFMlpANN.cxx:732
 MethodCFMlpANN.cxx:733
 MethodCFMlpANN.cxx:734
 MethodCFMlpANN.cxx:735
 MethodCFMlpANN.cxx:736
 MethodCFMlpANN.cxx:737
 MethodCFMlpANN.cxx:738
 MethodCFMlpANN.cxx:739
 MethodCFMlpANN.cxx:740
 MethodCFMlpANN.cxx:741
 MethodCFMlpANN.cxx:742
 MethodCFMlpANN.cxx:743
 MethodCFMlpANN.cxx:744
 MethodCFMlpANN.cxx:745
 MethodCFMlpANN.cxx:746
 MethodCFMlpANN.cxx:747
 MethodCFMlpANN.cxx:748
 MethodCFMlpANN.cxx:749
 MethodCFMlpANN.cxx:750
 MethodCFMlpANN.cxx:751
 MethodCFMlpANN.cxx:752
 MethodCFMlpANN.cxx:753
 MethodCFMlpANN.cxx:754