/*****************************************************************************

Copyright (c) 2005, Uwe Schmitt (http://www.procoders.net)
All rights reserved.

Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are met:

    * Redistributions of source code must retain the above copyright notice,
      this list of conditions and the following disclaimer.

    * Redistributions in binary form must reproduce the above copyright notice,
      this list of conditions and the following disclaimer in the documentation
      and/or other materials provided with the distribution.

    * Neither the name of procoders.net nor the names of its contributors may
      be used to endorse or promote products derived from this software without
      specific prior written permission.

THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR
ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
(INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

******************************************************************************/

#include "Classifier.hpp"
#include "QuadOptimData.hpp"
#include "Solver.hpp"

#include <iostream>
#include <stdexcept>
#include <cmath>
#include <limits>


using namespace pcsvm;

CSVMTrainInfo ClassifierCSVM::train(float cp, float cn)
{
        if (!hasKernel_) throw std::runtime_error("no kernel given.");
        if (!p_->size()) throw std::runtime_error("problem is empty.");
 
        unsigned int psize = p_->size();

        CSVMOptimData  da(p_, k_, cp, cn);

        for (unsigned int i=0; i<psize; ++i)
            y_[i]=da.y(i);

        // SVFloat alpha;

        //xx: todo: cache bei bedarf einbauen !
        //xx: todo: solver raus faktorisieren ! trotzdem standard löser anbieten !
        LibSVMSolver solver(&da, 1000000); 
        
        solver.setEps(1e-3);
        solver.doShrinking(true);
        solver.setMaxIter(100000);
        solver.solve(&alpha_); // schneller so, als alpha als return value

        //numSV_ = alpha_.entries();
        CSVMTrainInfo info;
        info.numSV = alpha_.entries();
        info.numIter = solver.solutionInfo.num_iter;

        for (SparseVector::const_iterator it = alpha_.begin(); it != alpha_.end(); ++it)
            info.svWeights[it->index]=it->value * y_[it->index];


        if (isLinear_)
        {
                // den einen support-vektor ausrechnen
                w_.clear();
                for (SparseVector::const_iterator it = alpha_.begin(); it != alpha_.end(); ++it)
                {
                    w_.plusAlphaYpsilon(it->value*y_[it->index], (*p_)[it->index].vector);
                        //w_ = w_ + (it->value) * y_[it->index]* ((p_->operator[])(it->index).vector);

                }
                //std::cout << "w: " <<w_ << std::endl;
        }
        else
        {
                //std::cout << "alpha: " << alpha_ << std::endl;
                supportVectors_.clear();
                // supportvectoren aus problem extrahieren
                for (SparseVector::const_iterator it = alpha_.begin(); it != alpha_.end(); ++it)
                        supportVectors_.push_back((*p_)[it->index].vector);
        }

        isTrained_ = true;

        // bias rho ausrechnen

        float sum=0;
        int   count = 0;

        float lb = std::numeric_limits<float>::min();
        float ub = std::numeric_limits<float>::max();

        for (SparseVector::const_iterator it= alpha_.begin(); it !=alpha_.end(); ++it)
        {
            int     i = it->index;
            float   alphai = it->value;
            //bool    yi = (*p_)[i].label;
            SparseVector xi = (*p_)[i].vector;

            float ci = y_[i] > 0 ? cp : cn;
            float yG = y_[i]-Phi(xi);

            if (alphai <= 0.0)
            {
                if (y_[i]>0)
                    ub= ub < yG ? ub : yG; 
                else
                    lb = lb > yG ? lb : yG;
            }
            else if (alphai >= ci)
            {
                if (y_[i]<0)
                    ub= ub < yG ? ub : yG; 
                else
                    lb = lb > yG ? lb : yG;
            }
            else// if (yi && 0.0f < alphai && alphai < cp || !yi && 0.0f < alphai && alphai < cn)
            {
                count += 1;
                sum += yG; // yi - Phi(xi);
            }
        }
        if (!count)
            // bin mit nicht sicher ob dieser fall ueberhaupt auftreten kann
            // throw std::runtime_error("severe implementation error: can not calculate bias rho due to missing support vectors. please contact author.");
            rho_ = (lb+ub)/2.0f;
        else
            rho_ = sum / count;

        // std::cout << rho_ << std::endl;
        
        return info;
}

CSVMTrainInfo ClassifierCSVM::train()
{
    float mval = std::numeric_limits<float>::max();
    return train(mval, mval);
}


CSVMTrainInfo ClassifierCSVM::train(float c)
{
        return train(c, c);
} 



float ClassifierCSVM::operator()(const SparseVector &what) const
{
        return Phi(what)+rho_;
}

float ClassifierCSVM::Phi(const SparseVector &what) const
// unbiased estimator, needed internally for computing bias !
{

        if (!hasKernel_) throw std::runtime_error("no kernel given.");
        if (!isTrained_) throw std::runtime_error("machine not trained yet");

        float val = 0.0;

        if (isLinear_) return (*k_)(what, w_);

        std::list<SparseVector>::const_iterator itSV  = supportVectors_.begin();
        SparseVector::const_iterator itAlpha = alpha_.begin();
        
        for (; itAlpha != alpha_.end(); ++itAlpha, ++itSV)
                val += y_[itAlpha->index] * itAlpha->value * (*k_)(what, *itSV);

        return val;
}

int ClassifierCSVM::getNumSVs() const
{
    if (!isTrained_) throw std::runtime_error("machine not trained yet");
    return alpha_.entries();
}

float ClassifierCSVM::getThreshold() const
{
    if (!isTrained_) throw std::runtime_error("machine not trained yet");
    return rho_;
}


//xx: TODO
//xx: TrainingINfo: #false positive, #false negatives in trainingsdaten + jeweiliger abstand !
//xx: abstand der supportvektoren im vergleich dazu  (2/||w||?)
/*
ClassifierCSVM::TrainingInfo ClassifierCSVM::calculateTrainingInfo()
{
        if (!hasKernel_) throw std::runtimeerror_("no kernel given.");
        if (!isTrained_) throw std::runtimeerror_("machine not trained yet");

        TrainingInfo info;

        info.numSV       = isLinear_ ? 1 : supportVectors_.size(); 
        info.problemSize = p_.size();
        info.numIter     = numIter_;

        const SVFloat& weightVec = isLinear_ ? alpha_ : alpha_;

        SVFloat::constiterator_ it;
        for (it = weightVec.begin(); it != weightVec.end(); ++it)
        {
            info.SVIndices.pushback_(it->index);
            info.SVWeights.pushback_(it->value);
        }

        //xx: folgendes ok ? gleicher wert fuer beide VZ ?
        float maxv = (float)-HUGEVAL_;
        float minv = (float)+HUGEVAL_;
        float tmp;

        for (ClassificationProblem::constiterator_ it = p_.begin(); it != p_.end(); ++it)
        {
                tmp = (*this)(it->vector);
                if (it->label && tmp<minv) 
                        minv=tmp;
                else if (!it->label && tmp>maxv)
                        maxv=tmp;
        }
        info.cp = minv;
        info.cn = -maxv; 

        return info;
}
*/