#ifndef __FILE_SNP_H_SEEN__
#define __FILE_SNP_H_SEEN__

/*-----------------------------------------------------------------------------

Copyright (C) 2004, 2006, 2009, 2014, 2016, 2019.

A. Ronald Gallant
Post Office Box 659
Chapel Hill NC 27514-0659
USA   

This program is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.

This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
GNU General Public License for more details.

You should have received a copy of the GNU General Public License along
with this program; if not, write to the Free Software Foundation, Inc.,
51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.

-----------------------------------------------------------------------------*/

#include "libsnp.h"
#include "snp_base.h"

namespace snp {

  class snpll;

  class control {
  private:
    std::string input_file;
    std::string output_file;
    REAL fnew;
    REAL fold;
    INTEGER nstart;
    INT_32BIT jseed;
    std::istream& control_stream;
  public:
    control(std::istream& ctrl_stream);
    bool read_line();
    std::string get_line();
    std::string get_input_file();
    std::string get_output_file();
    REAL get_fnew();
    REAL get_fold();
    INTEGER get_nstart();
    INT_32BIT get_jseed();
  };
    
  class parmfile {
  private:
    std::vector<std::string> history;
    optparms opm;
    datparms dpm;
    tranparms tpm;
    libsnp::snpden fn;
    libsnp::afunc afn;
    libsnp::ufunc ufn;
    libsnp::rfunc rfn;
    libsnp::afunc afn_mask;   // In these, -1 is fixed, and 0 is active.
    libsnp::ufunc ufn_mask;   // Add 1 to get the standard convention.
    libsnp::rfunc rfn_mask;
    bool rv; 
  public:
    parmfile();
    bool set_parms(const std::vector<std::string>& pfvec, std::ostream& detail);
    bool read_parms(const char* filename, std::ostream& detail);
    bool write_parms (const char* filename, std::string ctrl, snpll& ll) const;
    const optparms& get_optparms() const {return opm;}
    const datparms& get_datparms() const {return dpm;}
    const tranparms& get_tranparms() const {return tpm;}
    const libsnp::snpden& get_snpden() const {return fn;}
    const libsnp::afunc& get_afunc() const {return afn;}
    const libsnp::ufunc& get_ufunc() const {return ufn;}
    const libsnp::rfunc& get_rfunc() const {return rfn;}
    const libsnp::afunc& get_afunc_mask() const {return afn_mask;}
    const libsnp::ufunc& get_ufunc_mask() const {return ufn_mask;}
    const libsnp::rfunc& get_rfunc_mask() const {return rfn_mask;}
    void set_optparms(const optparms& op);
    void set_datparms(const datparms& dp);
    void set_tranparms(const tranparms& tp);
    void set_afunc(const libsnp::afunc& af);
    void set_ufunc(const libsnp::ufunc& uf);
    void set_rfunc(const libsnp::rfunc& rf);
  };
  
  class snpll {
  private:
    datparms dpm;
    const scl::realmat* Y;
    const scl::realmat* X;
    libsnp::snpden f;
    libsnp::afunc af;
    libsnp::ufunc uf;
    libsnp::rfunc rf;
    libsnp::afunc af_mask; //In these, -1 is fixed, and 0 is active.
    libsnp::ufunc uf_mask; //Add 1 to get the standard convention.
    libsnp::rfunc rf_mask;
    scl::realmat rho;
    scl::realmat theta;
    scl::intvec srvec;  //For MCMC retrict these theta elements to positive
    std::vector< std::pair<INTEGER,INTEGER> > rt;
    scl::realmat a0;
    scl::realmat A;
    scl::realmat b0;
    scl::realmat B;
    scl::realmat Rparms;
    bool compute_infmat;
    scl::realmat infmat;
    bool compute_scores;
    scl::realmat scores;
    void agglomerate();
    void distribute();
  public:
    snpll();
    snpll(datparms dp,
      libsnp::snpden fn, libsnp::afunc afn, 
      libsnp::ufunc ufn, libsnp::rfunc rfn,
      libsnp::afunc afn_mask, 
      libsnp::ufunc ufn_mask, libsnp::rfunc rfn_mask);
    snpll(const snpll& sl);
    ~snpll() { }
    snpll& operator=(const snpll& sl);
    void set_XY(const scl::realmat* x, const scl::realmat* y);  
    void set_rho(const scl::realmat& irho);
    void set_rho(const scl::realmat& irho, REAL fnew, 
      REAL fold, INT_32BIT& jseed);
    void set_n(INTEGER n) {dpm.n=n;}
    const scl::realmat* get_X() const { return X; }
    const scl::realmat* get_Y() const { return Y; }
    const libsnp::snpden& get_snpden() const {return f;}
    const libsnp::afunc& get_afunc() const {return af;}
    const libsnp::ufunc& get_ufunc() const {return uf;}
    const libsnp::rfunc& get_rfunc() const {return rf;}
    const datparms& get_datparms() const {return dpm;}
    const scl::realmat& get_rho() const {return rho;}
    const scl::realmat& get_theta() const {return theta;}
    const scl::intvec& get_srvec() const {return srvec;}
    REAL get_uf_stability() const {return uf.stability();}
    REAL get_rf_stability() const {return rf.stability();}
    const std::vector< std::pair<INTEGER,INTEGER> >& get_rt() const {return rt;}
    REAL log_likehood();
    REAL log_likehood(scl::realmat& dllwrho);
    REAL log_likehood(INTEGER& drop, scl::realmat& llv, scl::realmat& dllvwrho);
    REAL log_likehood(scl::realmat& dllwrho, scl::realmat& rho_infmat);
    REAL log_likehood(scl::realmat& dllwrho, scl::realmat& rho_infmat,
      scl::realmat& rho_scores);
    void rho_hessian(scl::realmat& rho_hessian);
    std::ostream& write_stats(std::ostream& os);
  };

  class sumry {
  private:
    std::ostream& os;
    REAL best_sn;
    std::string filename;
    std::string pname;
    const snp::snpll* ll_ptr;
    INTEGER n;
    REAL sn;
    INTEGER iter;
    INTEGER termination_code;
  public:
    sumry(std::ostream& sumry_os);
    void set_filename(const std::string& fn);
    void set_ll_ptr(const snp::snpll* llptr);
    void set_pname(const std::string ipname);
    void set_n(INTEGER in);
    void set_p(INTEGER ip); //Could get from ll, change later.
    void set_sn(REAL isn);  //Could get from ll but would destroy const.
    void set_iter(INTEGER itr);
    void set_termination_code(INTEGER code);
    bool write_header();
    bool write_line();
    bool write_partial_line();
  };
    
}

#endif
