source: Sophya/trunk/SophyaPI/PIext/nomtmatvecadapter.cc@ 3713

#include "sopnamsp.h"
#include "machdefs.h"
#include <stdlib.h>
#include <typeinfo>
#include <iostream>
#include <string>
#include <complex>

#include "datatype.h"

#include "tvector.h"
#include "objfitter.h"
#include "nomtmatvecadapter.h"
#include "piyfxdrw.h"
#include "pitvmaad.h"

#include "cimage.h"       // Pour la prise en charge des Image<T> (Avr 2007)
#include "pimgadapter.h"  // Adaptateur de tableaux PIImage pour Image<T>

#include "fioarr.h"

#include "nobjmgr.h"

//----------------------------------------------------------------
// Class Adaptateur d'objet (Pour NamedObjMgr) d'objet TMatrix<T> 
//----------------------------------------------------------------


/* --Methode-- */
template <class T> 
NOMAdapter_TMatrix<T>::NOMAdapter_TMatrix(TMatrix<T>* o)
  : NObjMgrAdapter(o)
{
mMtx = o;
}

/* --Methode-- */
template <class T> 
NOMAdapter_TMatrix<T>::~NOMAdapter_TMatrix()
{
}

/* --Methode-- */
template <class T> 
NObjMgrAdapter* NOMAdapter_TMatrix<T>::Clone(AnyDataObj* o)
{
TMatrix<T>* m = dynamic_cast<TMatrix<T> *>(o);
if (m) return ( new NOMAdapter_TMatrix<T>(m) );
return ( new NObjMgrAdapter(o) );
}

/* --Methode-- */
template <class T> 
string NOMAdapter_TMatrix<T>::GetDataObjType()
{
string type = "TMatrix< ";

TVector<T>* v = dynamic_cast<TVector<T> *>(mMtx);
if (v != NULL) type = "TVector< ";
else {
  Image<T>* img = dynamic_cast<Image<T> *>(mMtx);
  if (img != NULL) type = "Image< ";
}
// type +=  DecodeTypeIdName(typeid(T).name());
type += DataTypeInfo<T>::getTypeName();
type +=  " > ";
return(type);
}

/* --Methode-- */
template <class T> 
AnyDataObj* NOMAdapter_TMatrix<T>::CloneDataObj(bool share)
{
if (mMtx == NULL)  return(NULL);
TVector<T>* v = dynamic_cast<TVector<T> *>(mMtx);
if (v != NULL) return( new TVector<T>(*v, share) );
else {
  Image<T>* img = dynamic_cast<Image<T> *>(mMtx);
  if (img != NULL) return ( new Image<T>(*img, share) );
  else return ( new TMatrix<T>(*mMtx, share) );
}

}

/* --Methode-- */
template <class T> 
string NOMAdapter_TMatrix<T>::GetInfoString(vector<string>& opts)
{
  if (opts.size() == 0) return mMtx->InfoString();
  else {
    if (opts[0] == "rank") return string("2");
    else if (opts[0] == "sizes") {
      char buff[64];
      TVector<T>* v = dynamic_cast<TVector<T> *>(mMtx);
      if (v != NULL) sprintf(buff, "%ld", (long)mMtx->Size());
      else sprintf(buff, "%ld %ld", (long)mMtx->NRows(), (long)mMtx->NCols()); 
      return string(buff);
    }
    else if ((opts[0] == "size") || (opts[0] == "nelts")) {
      char buff[32];
      sprintf(buff, "%ld", (long)mMtx->Size()); 
      return string(buff);
    }
    else if ((opts[0] == "nrow") || (opts[0] == "nrows")) {
      char buff[32];
      sprintf(buff, "%ld", (long)mMtx->NRows()); 
      return string(buff);
    }
    else if ((opts[0] == "ncol") || (opts[0] == "ncols")) {
      char buff[32];
      sprintf(buff, "%ld", (long)mMtx->NCols()); 
      return string(buff);
    }
    else if (opts[0] == "sum") {
      MuTyV mtv(mMtx->Sum());
      string s;
      return mtv.Convert(s);
    }
    else if (opts[0] == "sumsq") {
      MuTyV mtv(mMtx->SumSq());
      string s;
      return mtv.Convert(s);
    }
    else if (opts[0] == "norm2") {
      MuTyV mtv(mMtx->Norm2());
      r_8 nrm2 = mtv.Convert(nrm2);  // we force the conversion to double
      mtv = nrm2;
      string s;
      return mtv.Convert(s);
    }
    else if ((opts[0] == "min")||(opts[0] == "max")||(opts[0] == "minmax")) {
      T amin, amax;
      MuTyV mtv;
      string s;
      mMtx->MinMax(amin, amax);
      if (opts[0] == "minmax") {
        mtv = amin;  mtv.Convert(s);
        s += " ";
        string mms;
        mtv = amax;  mtv.Convert(mms);
        s += mms;
      }
      else {
        if (opts[0] == "min")  mtv = amin;
        else if (opts[0] == "max")  mtv = amax;
        mtv.Convert(s); 
      }
      return s;
    }
    else if (opts[0] == "info") { // Acces aux valeurs stockes ds le DVList Info()
      if (opts.size() < 2)  return string(""); 
      else if (mMtx->Info().HasKey(opts[1]) == false) return string(""); 
      else return mMtx->Info().GetS(opts[1]);
    }
    else return "TMatrix.Att: rank size/nelts nrow/nrows ncol/ncols sum sumsq norm2 min max minmax info.varname";
  }
}

/* --Methode-- */
template <class T> 
int NOMAdapter_TMatrix<T>::PerformOperation(vector<string>& opts)
{
  bool ok = false;
  if (opts.size() >= 2) {
    int_4 kk = atoi(opts[1].c_str());
    TVector<T> * vrc = new TVector<T>;
    char buff[40];
    if ((opts[0] == "row") || (opts[0] == "line")) { 
      vrc->Share(mMtx->Row((sa_size_t)kk));  
      ok = true; 
      sprintf(buff,"row_%ld",(long)kk);
      cout << "PerformOperation(): Extracting Row(" << kk << ") from matrix" << endl;
    }
    else if ((opts[0] == "col") || (opts[0] == "column")) { 
      vrc->Share(mMtx->Column((sa_size_t)kk));  
      ok = true; 
      sprintf(buff,"col_%ld",(long)kk);
      cout << "PerformOperation(): Extracting Column(" << kk << ") from matrix" << endl;
    }
    if (ok) {
      string nvrc;
      if (opts.size() > 2)  nvrc = opts[2];
      else nvrc = buff;
      NamedObjMgr omg;
      omg.AddObj(vrc, nvrc, true); 
      return 0;
    }
  }
  
  cout << "NOMAdapter_TMatrix<T>::PerformOperation(): Error operation/arguments !" << endl; 
  cout << "...Valid args: row/line r [rowname] , col/column c [rowname]" << endl;
  return 1; 
}

/* --Methode-- */
template <class T> 
void NOMAdapter_TMatrix<T>::SavePPF(POutPersist& pos, string const & nom)
{
if (mMtx == NULL)  return;
TVector<T>* v = dynamic_cast<TVector<T> *>(mMtx);
if (v != NULL) {
  FIO_TArray<T> fio(v);
  fio.Write(pos, nom);
}
else {
  Image<T>* img = dynamic_cast<Image<T> *>(mMtx);
  if (img != NULL) {
    FIO_Image<T> fio(img);
    fio.Write(pos, nom);
  }
  else {
    FIO_TArray<T> fio(mMtx);
    fio.Write(pos, nom);
  }
}
return;
}

/* --Methode-- */
template <class T> 
void NOMAdapter_TMatrix<T>::Print(ostream& os, int lev)
{
if (lev < 3) mMtx->Show(os, false);
else  mMtx->Show(os, true);
if (lev > 0) mMtx->Print(os, 10*lev);
}

/* --Methode-- */
template <class T> 
PIDrawer* NOMAdapter_TMatrix<T>::GetDrawer(string & dopt)
{
TVector<T>* v = dynamic_cast<TVector<T> *>(mMtx);
if (v == NULL) return(NULL);
else {
  dopt = "thinline " + dopt;
  return( new PIYfXDrawer( new POTVectorAdapter<T>(v, false, StringToCmplxDispOption(dopt)), NULL, true) );
  }
}

/* --Methode-- */
template <class T> 
P2DArrayAdapter* NOMAdapter_TMatrix<T>::Get2DArray(string & dopt)
{
Image<T>* img = dynamic_cast<Image<T> *>(mMtx);
if (img != NULL) return ( new ImageAdapter<T>(img, false) );
 else return ( new POTMatrixAdapter<T>(mMtx, false, StringToCmplxDispOption(dopt)) );
}


// ---- Specialisation pour complexes -----
DECL_TEMP_SPEC  /* equivalent a template <> , pour SGI-CC en particulier */
P2DArrayAdapter* NOMAdapter_TMatrix< complex<r_4> >::Get2DArray(string & dopt)
{
  return ( new POTMatrixAdapter< complex<r_4> >(mMtx, false, StringToCmplxDispOption(dopt) ) );
}
DECL_TEMP_SPEC  /* equivalent a template <> , pour SGI-CC en particulier */
P2DArrayAdapter* NOMAdapter_TMatrix< complex<r_8> >::Get2DArray(string & dopt)
{
  return ( new POTMatrixAdapter< complex<r_8> >(mMtx, false, StringToCmplxDispOption(dopt) ) );
}
// -------------------------------------------------------------

/* --Methode-- */
template <class T> 
NTupleInterface* NOMAdapter_TMatrix<T>::GetNTupleInterface(bool& adel)
{
adel = true;
return( new NTupInt_TMatrix<T>(mMtx) );
}


/* --Methode-- */
template <class T> 
GeneralFitData* NOMAdapter_TMatrix<T>::GetGeneralFitData(bool& adel
 ,GeneralFitData::FitErrType errtype,double errscale,double errmin
 ,int i1,int i2,int j1,int j2)
{
adel = false;
if(!mMtx) return(NULL);

int nx,ny;
TVector<T>* vec = dynamic_cast<TVector<T> *>(mMtx);

if(vec) {  // C'est un TVector
  nx = vec->NElts();
  ny = 1;
} else {
  nx = mMtx->NRows();
  ny = mMtx->NCol();
  if(nx<=0 || ny<=0) return(NULL);
}

i1 = (i1<0||i1>=nx)? 0: i1;
i2 = (i2<0||i2>=nx||i2<i1)? nx-1: i2;
j1 = (j1<0||j1>=ny)? 0: j1;
j2 = (j2<0||j2>=ny||j2<j1)? ny-1: j2;

GeneralFitData* mGData = NULL;

if(vec) {  // C'est un TVector
  mGData = new GeneralFitData(1,(i2-i1+1),0);
  adel = true;
  for(int i=i1;i<=i2;i++) {
    double x = (double) i;
    double f = (*vec)(i);
    double e = 1.;
    e = GeneralFitData::ComputeError(f,e,errtype,errscale,errmin);
    mGData->AddData1(x,f,e);
  }
} else {
  mGData = new GeneralFitData(2,(i2-i1+1)*(j2-j1+1),0);
  adel = true;
  for(int i=i1;i<=i2;i++) for(int j=j1;j<=j2;j++) {
    double x = (double) i;
    double y = (double) j;
    double f = (*mMtx)(i,j);
    double e = 1.;
    e = GeneralFitData::ComputeError(f,e,errtype,errscale,errmin);
    mGData->AddData2(x,y,f,e);
  }
}

return mGData;
}

template <class T> 
AnyDataObj* NOMAdapter_TMatrix<T>::FitResidusObj(GeneralFit& mfit)
{
ArrayFitter<T> arrfit;
TVector<T>* vec = dynamic_cast<TVector<T> *>(mMtx);
if(vec) {
  TVector<T>* v = new TVector<T>(arrfit.FitResidus(*vec, mfit),true);
  return v;
} else {
  TMatrix<T>* m = new TMatrix<T>(arrfit.FitResidus(*mMtx,mfit),true);
  return m;
}
}

template <class T> 
AnyDataObj* NOMAdapter_TMatrix<T>::FitFunctionObj(GeneralFit& mfit)
{
ArrayFitter<T> arrfit;
TVector<T>* vec = dynamic_cast<TVector<T> *>(mMtx);
if(vec) {
  TVector<T>* v = new TVector<T>(arrfit.FitFunction(*vec,mfit),true);
  return v;
} else {
  TMatrix<T>* m = new TMatrix<T>(arrfit.FitFunction(*mMtx,mfit),true);
  return m;
}
}

// ---- Specialisation pour complexes -----
DECL_TEMP_SPEC  /* equivalent a template <> , pour SGI-CC en particulier */
GeneralFitData* NOMAdapter_TMatrix< complex<r_4> >::GetGeneralFitData(bool& adel
 ,GeneralFitData::FitErrType errtype,double errscale,double errmin
 ,int i1,int i2,int j1,int j2)
{
  return(NULL);
}

DECL_TEMP_SPEC  /* equivalent a template <> , pour SGI-CC en particulier */
AnyDataObj* NOMAdapter_TMatrix< complex<r_4> >::FitResidusObj(GeneralFit& mfit)
{
  return(NULL);
}

DECL_TEMP_SPEC  /* equivalent a template <> , pour SGI-CC en particulier */
AnyDataObj* NOMAdapter_TMatrix< complex<r_4> >::FitFunctionObj(GeneralFit& mfit)
{
  return(NULL);
}

DECL_TEMP_SPEC  /* equivalent a template <> , pour SGI-CC en particulier */
GeneralFitData* NOMAdapter_TMatrix< complex<r_8> >::GetGeneralFitData(bool& adel
 ,GeneralFitData::FitErrType errtype,double errscale,double errmin
 ,int i1,int i2,int j1,int j2)
{
  return(NULL);
}

DECL_TEMP_SPEC  /* equivalent a template <> , pour SGI-CC en particulier */
AnyDataObj* NOMAdapter_TMatrix< complex<r_8> >::FitResidusObj(GeneralFit& mfit)
{
  return(NULL);
}

DECL_TEMP_SPEC  /* equivalent a template <> , pour SGI-CC en particulier */
AnyDataObj* NOMAdapter_TMatrix< complex<r_8> >::FitFunctionObj(GeneralFit& mfit)
{
  return(NULL);
}

// -------------------------------------------------------------

/* --Methode-- */
template <class T> 
NTupInt_TMatrix<T>::NTupInt_TMatrix(TMatrix<T>* m)
{
mMtx = m;
}

/* --Methode-- */
template <class T> 
NTupInt_TMatrix<T>::~NTupInt_TMatrix()
{
}

/* --Methode-- */
template <class T> 
sa_size_t NTupInt_TMatrix<T>::NbLines() const
{
return( mMtx->NRows()*mMtx->NCols() );
}

/* --Methode-- */
template <class T> 
sa_size_t  NTupInt_TMatrix<T>::NbColumns() const 
{
return(8);
}

/* --Methode-- */
template <class T> 
r_8* NTupInt_TMatrix<T>::GetLineD(sa_size_t n) const 
{
int i,j;
if ((n < 0) || (n >= (int)(mMtx->NRows()*mMtx->NCols()) )) {
  mRet[0] = n;
  for(i=1; i<8; i++)  mRet[i] = 0.;
}
else { 
  i = n/mMtx->NCols(); j = n%mMtx->NCols(); 
  mRet[0] = n;  mRet[1] = i;  mRet[2] = j;  
  mRet[3] = (*mMtx)(i,j);
  mRet[4] = mRet[2];  mRet[5] = 0.;
  mRet[6] = mRet[2];  mRet[7] = 0.;
  }
return(mRet);
}

/* --Methode-- */
template <class T> 
string NTupInt_TMatrix<T>::VarList_C(const char* nx) const 
{
string nomx;
if (nx) nomx = nx;
else nomx = "_xh_";
string vardec = "double n,r,c,val,real,imag,mod,phas; \n";
vardec += "n = " + nomx + "[0];  r = " + nomx + "[1]; c = " + nomx + "[2]; \n"; 
vardec += "val = " + nomx + "[3]; \n";
vardec += "real = " + nomx + "[4];  imag = " + nomx + "[5]; \n";
vardec += "mod  = " + nomx + "[6];  phas = " + nomx + "[7]; \n";
return(vardec);
}

/* --Methode-- */
DECL_TEMP_SPEC  /* equivalent a template <> , pour SGI-CC en particulier */
r_8* NTupInt_TMatrix< complex<r_4> >::GetLineD(sa_size_t n) const 
{
int i,j;
if ((n < 0) || (n >= (int)(mMtx->NRows()*mMtx->NCols()) )) {
  mRet[0] = n;
  for(i=1; i<8; i++)  mRet[i] = 0.;
}
else { 
  i = n/mMtx->NCols(); j = n%mMtx->NCols(); 
  mRet[0] = n;   mRet[1] = i;  mRet[2] = j;  
  mRet[4] = (*mMtx)(i,j).real();  mRet[5] = (*mMtx)(i,j).imag();
  mRet[3] = mRet[6] = sqrt(mRet[4]*mRet[4]+mRet[5]*mRet[5]);
  mRet[7] = atan2(mRet[5], mRet[4]);
}
return(mRet);
}

DECL_TEMP_SPEC  /* equivalent a template <> , pour SGI-CC en particulier */
r_8* NTupInt_TMatrix< complex<r_8> >::GetLineD(sa_size_t n) const 
{
int i,j;
if ((n < 0) || (n >= (int)(mMtx->NRows()*mMtx->NCols()) )) {
  mRet[0] = n;
  for(i=1; i<8; i++)  mRet[i] = 0.;
}
else { 
  i = n/mMtx->NCols(); j = n%mMtx->NCols(); 
  mRet[0] = n;   mRet[1] = i;  mRet[2] = j;  
  mRet[4] = (*mMtx)(i,j).real();  mRet[5] = (*mMtx)(i,j).imag();
  mRet[3] = mRet[6] = sqrt(mRet[4]*mRet[4]+mRet[5]*mRet[5]);
  mRet[7] = atan2(mRet[5], mRet[4]);
}
return(mRet);
}


#ifdef __CXX_PRAGMA_TEMPLATES__
#pragma define_template NOMAdapter_TMatrix<uint_1>
#pragma define_template NOMAdapter_TMatrix<uint_2>
#pragma define_template NOMAdapter_TMatrix<uint_4>
#pragma define_template NOMAdapter_TMatrix<uint_8>
#pragma define_template NOMAdapter_TMatrix<int_1>
#pragma define_template NOMAdapter_TMatrix<int_2>
#pragma define_template NOMAdapter_TMatrix<int_4>
#pragma define_template NOMAdapter_TMatrix<int_8>
#pragma define_template NOMAdapter_TMatrix<r_4>
#pragma define_template NOMAdapter_TMatrix<r_8>
#pragma define_template NOMAdapter_TMatrix< complex<r_4> >
#pragma define_template NOMAdapter_TMatrix< complex<r_8> >
#pragma define_template NTupInt_TMatrix<uint_1>
#pragma define_template NTupInt_TMatrix<uint_2>
#pragma define_template NTupInt_TMatrix<uint_4>
#pragma define_template NTupInt_TMatrix<uint_8>
#pragma define_template NTupInt_TMatrix<int_1>
#pragma define_template NTupInt_TMatrix<int_2>
#pragma define_template NTupInt_TMatrix<int_4>
#pragma define_template NTupInt_TMatrix<int_8>
#pragma define_template NTupInt_TMatrix<r_4>
#pragma define_template NTupInt_TMatrix<r_8>
#pragma define_template NTupInt_TMatrix< complex<r_4> >
#pragma define_template NTupInt_TMatrix< complex<r_8> >
#endif

#if defined(ANSI_TEMPLATES) 
template class NOMAdapter_TMatrix<uint_1>;
template class NOMAdapter_TMatrix<uint_2>;
template class NOMAdapter_TMatrix<uint_4>;
template class NOMAdapter_TMatrix<uint_8>;
template class NOMAdapter_TMatrix<int_1>;
template class NOMAdapter_TMatrix<int_2>;
template class NOMAdapter_TMatrix<int_4>;
template class NOMAdapter_TMatrix<int_8>;
template class NOMAdapter_TMatrix<r_4>;
template class NOMAdapter_TMatrix<r_8>;
template class NOMAdapter_TMatrix< complex<r_4> >;
template class NOMAdapter_TMatrix< complex<r_8> >;
template class NTupInt_TMatrix<uint_1>;
template class NTupInt_TMatrix<uint_2>;
template class NTupInt_TMatrix<uint_4>;
template class NTupInt_TMatrix<uint_8>;
template class NTupInt_TMatrix<int_1>;
template class NTupInt_TMatrix<int_2>;
template class NTupInt_TMatrix<int_4>;
template class NTupInt_TMatrix<int_8>;
template class NTupInt_TMatrix<r_4>;
template class NTupInt_TMatrix<r_8>;
template class NTupInt_TMatrix< complex<r_4> >;
template class NTupInt_TMatrix< complex<r_8> >;
#endif