swh:1:snp:af87cd67498ef4fe47c76ed3e7caffe5b61facaf
Tip revision: 7dc0774baaef8b092ad0b1e3eaad7f482ad9a459 authored by Pere Mato on 23 April 2015, 16:30:25 UTC
Update ROOT version files to v5.34/30.
Update ROOT version files to v5.34/30.
Tip revision: 7dc0774
SVWorkingSet.cxx
// @(#)root/tmva $Id$
// Author: Andrzej Zemla
/**********************************************************************************
* Project: TMVA - a Root-integrated toolkit for multivariate data analysis *
* Package: TMVA *
* Class : SVWorkingSet *
* Web : http://tmva.sourceforge.net *
* *
* Description: *
* Implementation *
* *
* Authors (alphabetical): *
* Marcin Wolter <Marcin.Wolter@cern.ch> - IFJ PAN, Krakow, Poland *
* Andrzej Zemla <azemla@cern.ch> - IFJ PAN, Krakow, Poland *
* (IFJ PAN: Henryk Niewodniczanski Inst. Nucl. Physics, Krakow, Poland) *
* *
* Copyright (c) 2005: *
* CERN, Switzerland *
* MPI-K Heidelberg, Germany *
* PAN, Krakow, Poland *
* *
* 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) *
**********************************************************************************/
#include "TMath.h"
#include "TRandom3.h"
#ifndef ROOT_TMVA_MsgLogger
#include "TMVA/MsgLogger.h"
#endif
#include "TMVA/SVWorkingSet.h"
#include "TMVA/SVKernelFunction.h"
#include "TMVA/SVEvent.h"
#include "TMVA/SVKernelMatrix.h"
#include <vector>
#include <iostream>
//_______________________________________________________________________
TMVA::SVWorkingSet::SVWorkingSet()
: fdoRegression(kFALSE),
fInputData(0),
fSupVec(0),
fKFunction(0),
fKMatrix(0),
fTEventUp(0),
fTEventLow(0),
fB_low(1.),
fB_up(-1.),
fTolerance(0.01),
fLogger( new MsgLogger( "SVWorkingSet", kINFO ) )
{
// constructor
}
//_______________________________________________________________________
TMVA::SVWorkingSet::SVWorkingSet(std::vector<TMVA::SVEvent*>*inputVectors, SVKernelFunction* kernelFunction,
Float_t tol, Bool_t doreg)
: fdoRegression(doreg),
fInputData(inputVectors),
fSupVec(0),
fKFunction(kernelFunction),
fTEventUp(0),
fTEventLow(0),
fB_low(1.),
fB_up(-1.),
fTolerance(tol),
fLogger( new MsgLogger( "SVWorkingSet", kINFO ) )
{
// constructor
fKMatrix = new TMVA::SVKernelMatrix(inputVectors, kernelFunction);
Float_t *pt;
for( UInt_t i = 0; i < fInputData->size(); i++){
pt = fKMatrix->GetLine(i);
fInputData->at(i)->SetLine(pt);
fInputData->at(i)->SetNs(i);
if(fdoRegression) fInputData->at(i)->SetErrorCache(fInputData->at(i)->GetTarget());
}
TRandom3 rand;
UInt_t kk = rand.Integer(fInputData->size());
if(fdoRegression) {
fTEventLow = fTEventUp =fInputData->at(0);
fB_low = fTEventUp ->GetTarget() - fTolerance;
fB_up = fTEventLow->GetTarget() + fTolerance;
}
else{
while(1){
if(fInputData->at(kk)->GetTypeFlag()==-1){
fTEventLow = fInputData->at(kk);
break;
}
kk = rand.Integer(fInputData->size());
}
while (1){
if (fInputData->at(kk)->GetTypeFlag()==1) {
fTEventUp = fInputData->at(kk);
break;
}
kk = rand.Integer(fInputData->size());
}
}
fTEventUp ->SetErrorCache(fTEventUp->GetTarget());
fTEventLow->SetErrorCache(fTEventUp->GetTarget());
}
//_______________________________________________________________________
TMVA::SVWorkingSet::~SVWorkingSet()
{
// destructor
if (fKMatrix != 0) {delete fKMatrix; fKMatrix = 0;}
delete fLogger;
}
//_______________________________________________________________________
Bool_t TMVA::SVWorkingSet::ExamineExample( TMVA::SVEvent* jevt )
{
SVEvent* ievt=0;
Float_t fErrorC_J = 0.;
if( jevt->GetIdx()==0) fErrorC_J = jevt->GetErrorCache();
else{
Float_t *fKVals = jevt->GetLine();
fErrorC_J = 0.;
std::vector<TMVA::SVEvent*>::iterator idIter;
UInt_t k=0;
for(idIter = fInputData->begin(); idIter != fInputData->end(); idIter++){
if((*idIter)->GetAlpha()>0)
fErrorC_J += (*idIter)->GetAlpha()*(*idIter)->GetTypeFlag()*fKVals[k];
k++;
}
fErrorC_J -= jevt->GetTypeFlag();
jevt->SetErrorCache(fErrorC_J);
if((jevt->GetIdx() == 1) && (fErrorC_J < fB_up )){
fB_up = fErrorC_J;
fTEventUp = jevt;
}
else if ((jevt->GetIdx() == -1)&&(fErrorC_J > fB_low)) {
fB_low = fErrorC_J;
fTEventLow = jevt;
}
}
Bool_t converged = kTRUE;
if((jevt->GetIdx()>=0) && (fB_low - fErrorC_J > 2*fTolerance)) {
converged = kFALSE;
ievt = fTEventLow;
}
if((jevt->GetIdx()<=0) && (fErrorC_J - fB_up > 2*fTolerance)) {
converged = kFALSE;
ievt = fTEventUp;
}
if (converged) return kFALSE;
if(jevt->GetIdx()==0){
if(fB_low - fErrorC_J > fErrorC_J - fB_up) ievt = fTEventLow;
else ievt = fTEventUp;
}
if (TakeStep(ievt, jevt)) return kTRUE;
else return kFALSE;
}
//_______________________________________________________________________
Bool_t TMVA::SVWorkingSet::TakeStep(TMVA::SVEvent* ievt,TMVA::SVEvent* jevt )
{
if (ievt == jevt) return kFALSE;
std::vector<TMVA::SVEvent*>::iterator idIter;
const Float_t epsilon = 1e-8; //make it 1-e6 or 1-e5 to make it faster
Float_t type_I, type_J;
Float_t errorC_I, errorC_J;
Float_t alpha_I, alpha_J;
Float_t newAlpha_I, newAlpha_J;
Int_t s;
Float_t l, h, lobj = 0, hobj = 0;
Float_t eta;
type_I = ievt->GetTypeFlag();
alpha_I = ievt->GetAlpha();
errorC_I = ievt->GetErrorCache();
type_J = jevt->GetTypeFlag();
alpha_J = jevt->GetAlpha();
errorC_J = jevt->GetErrorCache();
s = Int_t( type_I * type_J );
Float_t c_i = ievt->GetCweight();
Float_t c_j = jevt->GetCweight();
// compute l, h objective function
if (type_I == type_J) {
Float_t gamma = alpha_I + alpha_J;
if ( c_i > c_j ) {
if ( gamma < c_j ) {
l = 0;
h = gamma;
}
else{
h = c_j;
if ( gamma < c_i )
l = 0;
else
l = gamma - c_i;
}
}
else {
if ( gamma < c_i ){
l = 0;
h = gamma;
}
else {
l = gamma - c_i;
if ( gamma < c_j )
h = gamma;
else
h = c_j;
}
}
}
else {
Float_t gamma = alpha_I - alpha_J;
if (gamma > 0) {
l = 0;
if ( gamma >= (c_i - c_j) )
h = c_i - gamma;
else
h = c_j;
}
else {
l = -gamma;
if ( (c_i - c_j) >= gamma)
h = c_j;
else
h = c_i - gamma;
}
}
if (l == h) return kFALSE;
Float_t kernel_II, kernel_IJ, kernel_JJ;
kernel_II = fKMatrix->GetElement(ievt->GetNs(),ievt->GetNs());
kernel_IJ = fKMatrix->GetElement(ievt->GetNs(), jevt->GetNs());
kernel_JJ = fKMatrix->GetElement(jevt->GetNs(),jevt->GetNs());
eta = 2*kernel_IJ - kernel_II - kernel_JJ;
if (eta < 0) {
newAlpha_J = alpha_J + (type_J*( errorC_J - errorC_I ))/eta;
if (newAlpha_J < l) newAlpha_J = l;
else if (newAlpha_J > h) newAlpha_J = h;
}
else {
Float_t c_I = eta/2;
Float_t c_J = type_J*( errorC_I - errorC_J ) - eta * alpha_J;
lobj = c_I * l * l + c_J * l;
hobj = c_I * h * h + c_J * h;
if (lobj > hobj + epsilon) newAlpha_J = l;
else if (lobj < hobj - epsilon) newAlpha_J = h;
else newAlpha_J = alpha_J;
}
if (TMath::Abs( newAlpha_J - alpha_J ) < ( epsilon * ( newAlpha_J + alpha_J+ epsilon ))){
return kFALSE;
//it spends here to much time... it is stupido
}
newAlpha_I = alpha_I - s*( newAlpha_J - alpha_J );
if (newAlpha_I < 0) {
newAlpha_J += s* newAlpha_I;
newAlpha_I = 0;
}
else if (newAlpha_I > c_i) {
Float_t temp = newAlpha_I - c_i;
newAlpha_J += s * temp;
newAlpha_I = c_i;
}
Float_t dL_I = type_I * ( newAlpha_I - alpha_I );
Float_t dL_J = type_J * ( newAlpha_J - alpha_J );
Int_t k = 0;
for(idIter = fInputData->begin(); idIter != fInputData->end(); idIter++){
k++;
if((*idIter)->GetIdx()==0){
Float_t ii = fKMatrix->GetElement(ievt->GetNs(), (*idIter)->GetNs());
Float_t jj = fKMatrix->GetElement(jevt->GetNs(), (*idIter)->GetNs());
(*idIter)->UpdateErrorCache(dL_I * ii + dL_J * jj);
}
}
ievt->SetAlpha(newAlpha_I);
jevt->SetAlpha(newAlpha_J);
// set new indexes
SetIndex(ievt);
SetIndex(jevt);
// update error cache
ievt->SetErrorCache(errorC_I + dL_I*kernel_II + dL_J*kernel_IJ);
jevt->SetErrorCache(errorC_J + dL_I*kernel_IJ + dL_J*kernel_JJ);
// compute fI_low, fB_low
fB_low = -1*1e30;
fB_up = 1e30;
for(idIter = fInputData->begin(); idIter != fInputData->end(); idIter++){
if((*idIter)->GetIdx()==0){
if((*idIter)->GetErrorCache()> fB_low){
fB_low = (*idIter)->GetErrorCache();
fTEventLow = (*idIter);
}
if( (*idIter)->GetErrorCache()< fB_up){
fB_up =(*idIter)->GetErrorCache();
fTEventUp = (*idIter);
}
}
}
// for optimized alfa's
if (fB_low < TMath::Max(ievt->GetErrorCache(), jevt->GetErrorCache())) {
if (ievt->GetErrorCache() > fB_low) {
fB_low = ievt->GetErrorCache();
fTEventLow = ievt;
}
else {
fB_low = jevt->GetErrorCache();
fTEventLow = jevt;
}
}
if (fB_up > TMath::Max(ievt->GetErrorCache(), jevt->GetErrorCache())) {
if (ievt->GetErrorCache()< fB_low) {
fB_up =ievt->GetErrorCache();
fTEventUp = ievt;
}
else {
fB_up =jevt->GetErrorCache() ;
fTEventUp = jevt;
}
}
return kTRUE;
}
//_______________________________________________________________________
Bool_t TMVA::SVWorkingSet::Terminated()
{
if((fB_up > fB_low - 2*fTolerance)) return kTRUE;
return kFALSE;
}
//_______________________________________________________________________
void TMVA::SVWorkingSet::Train(UInt_t nMaxIter)
{
// train the SVM
Int_t numChanged = 0;
Int_t examineAll = 1;
Float_t numChangedOld = 0;
Int_t deltaChanges = 0;
UInt_t numit = 0;
std::vector<TMVA::SVEvent*>::iterator idIter;
while ((numChanged > 0) || (examineAll > 0)) {
numChanged = 0;
if (examineAll) {
for (idIter = fInputData->begin(); idIter!=fInputData->end(); idIter++){
if(!fdoRegression) numChanged += (UInt_t)ExamineExample(*idIter);
else numChanged += (UInt_t)ExamineExampleReg(*idIter);
}
}
else {
for (idIter = fInputData->begin(); idIter!=fInputData->end(); idIter++) {
if ((*idIter)->IsInI0()) {
if(!fdoRegression) numChanged += (UInt_t)ExamineExample(*idIter);
else numChanged += (UInt_t)ExamineExampleReg(*idIter);
if (Terminated()) {
numChanged = 0;
break;
}
}
}
}
if (examineAll == 1) examineAll = 0;
else if (numChanged == 0 || numChanged < 10 || deltaChanges > 3 ) examineAll = 1;
if (numChanged == numChangedOld) deltaChanges++;
else deltaChanges = 0;
numChangedOld = numChanged;
++numit;
if (numit >= nMaxIter) {
*fLogger << kWARNING
<< "Max number of iterations exceeded. "
<< "Training may not be completed. Try use less Cost parameter" << Endl;
break;
}
}
}
//_______________________________________________________________________
void TMVA::SVWorkingSet::SetIndex( TMVA::SVEvent* event )
{
if( (0< event->GetAlpha()) && (event->GetAlpha()< event->GetCweight()))
event->SetIdx(0);
if( event->GetTypeFlag() == 1){
if( event->GetAlpha() == 0)
event->SetIdx(1);
else if( event->GetAlpha() == event->GetCweight() )
event->SetIdx(-1);
}
if( event->GetTypeFlag() == -1){
if( event->GetAlpha() == 0)
event->SetIdx(-1);
else if( event->GetAlpha() == event->GetCweight() )
event->SetIdx(1);
}
}
//_______________________________________________________________________
void TMVA::SVWorkingSet::PrintStat()
{
std::vector<TMVA::SVEvent*>::iterator idIter;
UInt_t counter = 0;
for( idIter = fInputData->begin(); idIter != fInputData->end(); idIter++)
if((*idIter)->GetAlpha() !=0) counter++;
}
//_______________________________________________________________________
std::vector<TMVA::SVEvent*>* TMVA::SVWorkingSet::GetSupportVectors()
{
std::vector<TMVA::SVEvent*>::iterator idIter;
if( fSupVec != 0) {delete fSupVec; fSupVec = 0; }
fSupVec = new std::vector<TMVA::SVEvent*>(0);
for( idIter = fInputData->begin(); idIter != fInputData->end(); idIter++){
if((*idIter)->GetDeltaAlpha() !=0){
fSupVec->push_back((*idIter));
}
}
return fSupVec;
}
//for regression
Bool_t TMVA::SVWorkingSet::TakeStepReg(TMVA::SVEvent* ievt,TMVA::SVEvent* jevt )
{
if (ievt == jevt) return kFALSE;
std::vector<TMVA::SVEvent*>::iterator idIter;
const Float_t epsilon = 0.001*fTolerance;//TODO
const Float_t kernel_II = fKMatrix->GetElement(ievt->GetNs(),ievt->GetNs());
const Float_t kernel_IJ = fKMatrix->GetElement(ievt->GetNs(),jevt->GetNs());
const Float_t kernel_JJ = fKMatrix->GetElement(jevt->GetNs(),jevt->GetNs());
//compute eta & gamma
const Float_t eta = -2*kernel_IJ + kernel_II + kernel_JJ;
const Float_t gamma = ievt->GetDeltaAlpha() + jevt->GetDeltaAlpha();
//TODO CHECK WHAT IF ETA <0
//w.r.t Mercer's conditions it should never happen, but what if?
Bool_t caseA, caseB, caseC, caseD, terminated;
caseA = caseB = caseC = caseD = terminated = kFALSE;
Float_t b_alpha_i, b_alpha_j, b_alpha_i_p, b_alpha_j_p; //temporary lagrange multipliers
const Float_t b_cost_i = ievt->GetCweight();
const Float_t b_cost_j = jevt->GetCweight();
b_alpha_i = ievt->GetAlpha();
b_alpha_j = jevt->GetAlpha();
b_alpha_i_p = ievt->GetAlpha_p();
b_alpha_j_p = jevt->GetAlpha_p();
//calculate deltafi
Float_t deltafi = ievt->GetErrorCache()-jevt->GetErrorCache();
// main loop
while(!terminated) {
const Float_t null = 0.; //!!! dummy float null declaration because of problems with TMath::Max/Min(Float_t, Float_t) function
Float_t low, high;
Float_t tmp_alpha_i, tmp_alpha_j;
tmp_alpha_i = tmp_alpha_j = 0.;
//TODO check this conditions, are they proper
if((caseA == kFALSE) && (b_alpha_i > 0 || (b_alpha_i_p == 0 && deltafi > 0)) && (b_alpha_j > 0 || (b_alpha_j_p == 0 && deltafi < 0)))
{
//compute low, high w.r.t a_i, a_j
low = TMath::Max( null, gamma - b_cost_j );
high = TMath::Min( b_cost_i , gamma);
if(low<high){
tmp_alpha_j = b_alpha_j - (deltafi/eta);
tmp_alpha_j = TMath::Min(tmp_alpha_j,high );
tmp_alpha_j = TMath::Max(low ,tmp_alpha_j);
tmp_alpha_i = b_alpha_i - (tmp_alpha_j - b_alpha_j);
//update Li & Lj if change is significant (??)
if( IsDiffSignificant(b_alpha_j,tmp_alpha_j, epsilon) || IsDiffSignificant(b_alpha_i,tmp_alpha_i, epsilon)){
b_alpha_j = tmp_alpha_j;
b_alpha_i = tmp_alpha_i;
}
}
else
terminated = kTRUE;
caseA = kTRUE;
}
else if((caseB==kFALSE) && (b_alpha_i>0 || (b_alpha_i_p==0 && deltafi >2*epsilon )) && (b_alpha_j_p>0 || (b_alpha_j==0 && deltafi>2*epsilon)))
{
//compute LH w.r.t. a_i, a_j*
low = TMath::Max( null, gamma ); //TODO
high = TMath::Min( b_cost_i , b_cost_j + gamma);
if(low<high){
tmp_alpha_j = b_alpha_j_p - ((deltafi-2*epsilon)/eta);
tmp_alpha_j = TMath::Min(tmp_alpha_j,high);
tmp_alpha_j = TMath::Max(low,tmp_alpha_j);
tmp_alpha_i = b_alpha_i - (tmp_alpha_j - b_alpha_j_p);
//update alphai alphaj_p
if( IsDiffSignificant(b_alpha_j_p,tmp_alpha_j, epsilon) || IsDiffSignificant(b_alpha_i,tmp_alpha_i, epsilon)){
b_alpha_j_p = tmp_alpha_j;
b_alpha_i = tmp_alpha_i;
}
}
else
terminated = kTRUE;
caseB = kTRUE;
}
else if((caseC==kFALSE) && (b_alpha_i_p>0 || (b_alpha_i==0 && deltafi < -2*epsilon )) && (b_alpha_j>0 || (b_alpha_j_p==0 && deltafi< -2*epsilon)))
{
//compute LH w.r.t. alphai_p alphaj
low = TMath::Max(null, -gamma );
high = TMath::Min(b_cost_i, -gamma+b_cost_j);
if(low<high){
tmp_alpha_j = b_alpha_j - ((deltafi+2*epsilon)/eta);
tmp_alpha_j = TMath::Min(tmp_alpha_j,high );
tmp_alpha_j = TMath::Max(low ,tmp_alpha_j);
tmp_alpha_i = b_alpha_i_p - (tmp_alpha_j - b_alpha_j);
//update alphai_p alphaj
if( IsDiffSignificant(b_alpha_j,tmp_alpha_j, epsilon) || IsDiffSignificant(b_alpha_i_p,tmp_alpha_i, epsilon)){
b_alpha_j = tmp_alpha_j;
b_alpha_i_p = tmp_alpha_i;
}
}
else
terminated = kTRUE;
caseC = kTRUE;
}
else if((caseD == kFALSE) &&
(b_alpha_i_p>0 || (b_alpha_i==0 && deltafi <0 )) &&
(b_alpha_j_p>0 || (b_alpha_j==0 && deltafi >0 )))
{
//compute LH w.r.t. alphai_p alphaj_p
low = TMath::Max(null,-gamma - b_cost_j);
high = TMath::Min(b_cost_i, -gamma);
if(low<high){
tmp_alpha_j = b_alpha_j_p + (deltafi/eta);
tmp_alpha_j = TMath::Min(tmp_alpha_j,high );
tmp_alpha_j = TMath::Max(low ,tmp_alpha_j);
tmp_alpha_i = b_alpha_i_p - (tmp_alpha_j - b_alpha_j_p);
if( IsDiffSignificant(b_alpha_j_p,tmp_alpha_j, epsilon) || IsDiffSignificant(b_alpha_i_p,tmp_alpha_i, epsilon)){
b_alpha_j_p = tmp_alpha_j;
b_alpha_i_p = tmp_alpha_i;
}
}
else
terminated = kTRUE;
caseD = kTRUE;
}
else
terminated = kTRUE;
}
// TODO ad commment how it was calculated
deltafi += ievt->GetDeltaAlpha()*(kernel_II - kernel_IJ) + jevt->GetDeltaAlpha()*(kernel_IJ - kernel_JJ);
if( IsDiffSignificant(b_alpha_i, ievt->GetAlpha(), epsilon) ||
IsDiffSignificant(b_alpha_j, jevt->GetAlpha(), epsilon) ||
IsDiffSignificant(b_alpha_i_p, ievt->GetAlpha_p(), epsilon) ||
IsDiffSignificant(b_alpha_j_p, jevt->GetAlpha_p(), epsilon) ){
//TODO check if these conditions might be easier
//TODO write documentation for this
const Float_t diff_alpha_i = ievt->GetDeltaAlpha()+b_alpha_i_p - ievt->GetAlpha();
const Float_t diff_alpha_j = jevt->GetDeltaAlpha()+b_alpha_j_p - jevt->GetAlpha();
//update error cache
Int_t k = 0;
for(idIter = fInputData->begin(); idIter != fInputData->end(); idIter++){
k++;
//there will be some changes in Idx notation
if((*idIter)->GetIdx()==0){
Float_t k_ii = fKMatrix->GetElement(ievt->GetNs(), (*idIter)->GetNs());
Float_t k_jj = fKMatrix->GetElement(jevt->GetNs(), (*idIter)->GetNs());
(*idIter)->UpdateErrorCache(diff_alpha_i * k_ii + diff_alpha_j * k_jj);
}
}
//store new alphas in SVevents
ievt->SetAlpha(b_alpha_i);
jevt->SetAlpha(b_alpha_j);
ievt->SetAlpha_p(b_alpha_i_p);
jevt->SetAlpha_p(b_alpha_j_p);
//TODO update Idexes
// compute fI_low, fB_low
fB_low = -1*1e30;
fB_up =1e30;
for(idIter = fInputData->begin(); idIter != fInputData->end(); idIter++){
if((!(*idIter)->IsInI3()) && ((*idIter)->GetErrorCache()> fB_low)){
fB_low = (*idIter)->GetErrorCache();
fTEventLow = (*idIter);
}
if((!(*idIter)->IsInI2()) && ((*idIter)->GetErrorCache()< fB_up)){
fB_up =(*idIter)->GetErrorCache();
fTEventUp = (*idIter);
}
}
return kTRUE;
} else return kFALSE;
}
//_______________________________________________________________________
Bool_t TMVA::SVWorkingSet::ExamineExampleReg(TMVA::SVEvent* jevt)
{
Float_t feps = 1e-7;// TODO check which value is the best
SVEvent* ievt=0;
Float_t fErrorC_J = 0.;
if( jevt->IsInI0()) {
fErrorC_J = jevt->GetErrorCache();
}
else{
Float_t *fKVals = jevt->GetLine();
fErrorC_J = 0.;
std::vector<TMVA::SVEvent*>::iterator idIter;
UInt_t k=0;
for(idIter = fInputData->begin(); idIter != fInputData->end(); idIter++){
fErrorC_J -= (*idIter)->GetDeltaAlpha()*fKVals[k];
k++;
}
fErrorC_J += jevt->GetTarget();
jevt->SetErrorCache(fErrorC_J);
if(jevt->IsInI1()){
if(fErrorC_J + feps < fB_up ){
fB_up = fErrorC_J + feps;
fTEventUp = jevt;
}
else if(fErrorC_J -feps > fB_low) {
fB_low = fErrorC_J - feps;
fTEventLow = jevt;
}
}else if((jevt->IsInI2()) && (fErrorC_J + feps > fB_low)){
fB_low = fErrorC_J + feps;
fTEventLow = jevt;
}else if((jevt->IsInI3()) && (fErrorC_J - feps < fB_up)){
fB_up = fErrorC_J - feps;
fTEventUp = jevt;
}
}
Bool_t converged = kTRUE;
//case 1
if(jevt->IsInI0a()){
if( fB_low -fErrorC_J + feps > 2*fTolerance){
converged = kFALSE;
ievt = fTEventLow;
if(fErrorC_J-feps-fB_up > fB_low-fErrorC_J+feps){
ievt = fTEventUp;
}
}else if(fErrorC_J -feps - fB_up > 2*fTolerance){
converged = kFALSE;
ievt = fTEventUp;
if(fB_low - fErrorC_J+feps > fErrorC_J-feps -fB_up){
ievt = fTEventLow;
}
}
}
//case 2
if(jevt->IsInI0b()){
if( fB_low -fErrorC_J - feps > 2*fTolerance){
converged = kFALSE;
ievt = fTEventLow;
if(fErrorC_J+feps-fB_up > fB_low-fErrorC_J-feps){
ievt = fTEventUp;
}
}else if(fErrorC_J + feps - fB_up > 2*fTolerance){
converged = kFALSE;
ievt = fTEventUp;
if(fB_low - fErrorC_J-feps > fErrorC_J+feps -fB_up){
ievt = fTEventLow;
}
}
}
//case 3
if(jevt->IsInI1()){
if( fB_low -fErrorC_J - feps > 2*fTolerance){
converged = kFALSE;
ievt = fTEventLow;
if(fErrorC_J+feps-fB_up > fB_low-fErrorC_J-feps){
ievt = fTEventUp;
}
}else if(fErrorC_J - feps - fB_up > 2*fTolerance){
converged = kFALSE;
ievt = fTEventUp;
if(fB_low - fErrorC_J+feps > fErrorC_J-feps -fB_up){
ievt = fTEventLow;
}
}
}
//case 4
if(jevt->IsInI2()){
if( fErrorC_J + feps -fB_up > 2*fTolerance){
converged = kFALSE;
ievt = fTEventUp;
}
}
//case 5
if(jevt->IsInI3()){
if(fB_low -fErrorC_J +feps > 2*fTolerance){
converged = kFALSE;
ievt = fTEventLow;
}
}
if(converged) return kFALSE;
if (TakeStepReg(ievt, jevt)) return kTRUE;
else return kFALSE;
}
Bool_t TMVA::SVWorkingSet::IsDiffSignificant(Float_t a_i, Float_t a_j, Float_t eps)
{
if( TMath::Abs(a_i - a_j) > eps*(a_i + a_j + eps)) return kTRUE;
else return kFALSE;
}