/* MELODIC - Multivariate exploratory linear optimized decomposition into independent components melgmix.h - class for Gaussian/Gamma Mixture Model Christian F. Beckmann, FMRIB Analysis Group Copyright (C) 1999-2013 University of Oxford */ /* CCOPYRIGHT */ #ifndef __MELGMIX_h #define __MELGMIX_h #include "newimage/newimageall.h" #include "utils/log.h" #include "melodic.h" #include "utils/options.h" #include "meloptions.h" //#include "melreport.h" using namespace Utilities; using namespace NEWIMAGE; namespace Melodic{ class MelGMix{ public: MelGMix(MelodicOptions &popts, Log &plogger): opts(popts), logger(plogger){} ~MelGMix() { //mainhtml << endl << "<hr></CENTER></BODY></HTML>" << endl; } void save(); void setup(const RowVector& dat, const string dirname, int here, volume<float> themask, volume<float> themean, int num_mix = 3, float eps = 0.0, bool fixdim = false); void gmmfit(); void ggmfit(); inline void fit(string mtype = string("GGM")){ mmtype = mtype; if(mmtype==string("GGM")) this->ggmfit(); else this->gmmfit(); //re-insert mean and stdev data = data*datastdev + datamean; //threshmaps = threshmaps*datastdev + datamean; means = means*datastdev + datamean; vars = vars*datastdev*datastdev; } inline Matrix threshold(string levels){ return this->threshold(data, levels); } inline Matrix threshold(RowVector& levels){ return this->threshold(data, levels); } Matrix threshold(const RowVector& dat, Matrix& levels); Matrix threshold(const RowVector& dat, string levels); void status(const string &txt); inline RowVector& get_means() {return means;} inline void set_means(RowVector& Arg) {means = Arg;} inline RowVector& get_vars() {return vars;} inline void set_vars(RowVector& Arg) {vars = Arg;} inline RowVector& get_pi() {return props;} inline void set_pi(RowVector& Arg) {props = Arg;} inline RowVector& get_data() {return data;} inline void set_data(RowVector& Arg) {data = Arg;} inline RowVector& get_prob() {return probmap;} inline float get_eps() {return epsilon;} inline void set_eps(float Arg) {epsilon = Arg;} inline Matrix& get_threshmaps() {return threshmaps;} inline void set_threshmaps(Matrix& Arg) {threshmaps = Arg;} inline bool isfitted(){return fitted;} inline int mixtures(){return nummix;} inline string get_type() { return mmtype;} inline void set_type(string Arg) { mmtype = Arg;} inline string get_prefix() { return prefix;} inline void set_prefix(string Arg) { prefix = Arg;} inline RowVector get_probmap() {return probmap;} inline float get_offset() {return offset;} inline void set_offset(float Arg) {offset = Arg;} inline void flipres(int num){ means = -means; data = -data; threshmaps = -threshmaps; if(mmtype=="GGM"){ float tmp; tmp= means(2);means(2)=means(3);means(3)=tmp; tmp=vars(2);vars(2)=vars(3);vars(3)=tmp; tmp=props(2);props(2)=props(3);props(3)=tmp; } } void create_rep(); inline void add_infstr(string what){ threshinfo.push_back(what); } inline string get_infstr(int num){ if((threshinfo.size()<(unsigned int)(num-1))||(num<1)) return string(""); else return threshinfo[num-1]; } inline int size_infstr(){ return threshinfo.size(); } inline void clear_infstr(){ threshinfo.clear(); } inline void smooth_probs(float howmuch){ volume4D<float> tempVol; tempVol.setmatrix(probmap,Mask); tempVol[0]= smooth(tempVol[0],howmuch); probmap = tempVol.matrix(Mask); } inline Matrix get_params(){ Matrix tmp = zeros(3,means.Ncols()); tmp.Row(1) = means; tmp.Row(2) = vars; tmp.Row(3) = props; return tmp; } double datamean; double datastdev; private: __attribute__((unused)) MelodicOptions &opts; __attribute__((unused)) Log &logger; //global log file //Log mainhtml; void gmmupdate(); float gmmevidence(); void gmmreducemm(); void add_params(Matrix& mu, Matrix& sig, Matrix& pi, float logLH, float MDL, float Evi, bool advance = false); void get_params(int index, Matrix& mu, Matrix& sig, Matrix& pi, float logLH, float MDL, float Evi); Matrix Params; Matrix threshmaps; RowVector means; RowVector vars; RowVector props; RowVector data; RowVector probmap; volume<float> Mean; volume<float> Mask; float epsilon; float logprobY; float MDL; float Evi; float offset; int nummix; int numdata; int cnumber; bool fitted; bool fixdim; string prefix; string mmtype; string dirname; vector<string> threshinfo; }; } #endif