source: Sophya/trunk/SophyaLib/NTools/generaldata.h@ 3853

// This may look like C code, but it is really -*- C++ -*-
//                         C.Magneville          04/99
#ifndef GENERALDATA_SEEN
#define GENERALDATA_SEEN

#include "objfio.h"
#include <iostream>
#include "pexceptions.h"
#include "ppersist.h"
#include "ntupintf.h"
#include "poly.h"

namespace SOPHYA {

//  Forward class declaration for Fits handler
template <class T>  class FitsHandler;

class GeneralFit;

//================================================================
// GENERALFITDATA
//================================================================

//! Classe de stoquage de donnees a plusieurs variables avec erreur
class GeneralFitData : public AnyDataObj , public NTupleInterface {
  friend class GeneralFit;
  friend class ObjFileIO<GeneralFitData>;
  friend class FitsHandler<GeneralFitData>;
public:
  //! Valeurs par defaut pour l'utilisation des erreurs
  enum {
    Def_ErrF = 1, //!< erreurs sur F par defaut
    Def_ErrX = 0  //!< pas d'erreurs sur les variables X,Y,... par defaut
    };

  GeneralFitData(uint_4 nVar, uint_4 nDataAlloc, uint_2 errx=0);
  GeneralFitData(const GeneralFitData& data, bool clean=false);
  GeneralFitData();
  virtual ~GeneralFitData();

  void SetDataPtr(int_4 ptr = 0);

  void KillData(int_4 i);
  void KillData(int_4 i1,int_4 i2);

  void ValidData(int_4 i);
  void ValidData(int_4 i1,int_4 i2);
  void ValidData();

  void RedefineData1(int_4 i,double x,double f,double err=Def_ErrF,double errx=Def_ErrX);
  void RedefineData2(int_4 i,double x,double y,double f,double err=Def_ErrF
                    ,double errx=Def_ErrX,double erry=Def_ErrX);
  void RedefineData(int_4 i,double* xp,double f,double err=Def_ErrF,double* errxp=NULL);

  void AddData1(double x, double f, double err=Def_ErrF,double errx=Def_ErrX);
  void AddData2(double x, double y, double f, double err=Def_ErrF
               ,double errx=Def_ErrX,double erry=Def_ErrX);
  void AddData(double* xp, double f, double err=Def_ErrF,double* errxp=NULL);
  void AddData(float* xp, float f, float err=Def_ErrF,float* errxp=NULL);

  void SetData1(int_4 nData,double* x,double* f,double *err=NULL,double *errx=NULL);
  void SetData1(int_4 nData,float* x,float* f,float* err=NULL,float *errx=NULL);
  void SetData2(int_4 nData,double* x,double* y,double* f,double *err=NULL
               ,double *errx=NULL,double *erry=NULL);
  void SetData2(int_4 nData,float* x,float* y,float* f,float* err=NULL
               ,float *errx=NULL,float *erry=NULL);
  void SetData(int_4 nData,double** xp,double *f,double *err=NULL,double** errxp=NULL);
  void SetData(int_4 nData,float** xp,float* f,float* err=NULL,float** errxp=NULL);

  void PrintStatus() const;
  void PrintData(int_4 i) const;
  void PrintData(int_4 i1,int_4 i2) const;
  void PrintData() const;
  void Show(ostream& os) const;
  inline void Show() const {Show(cout);}

  //! Operator equal
  GeneralFitData& operator = (const GeneralFitData& g);

  //! Retourne la place restante dans la structure (nombre de donnees que l'on peut encore stoquer).
  inline int_4 GetSpaceFree() const { return mNDataAlloc - mNData; }
  //! Retourne le nombre de variables Xi
  inline int_4 NVar()       const {return mNVar;}
  //! Retourne le nombre de donnees
  inline int_4 NData()      const {return mNData;}
  //! Retourne le nombre de bonnes donnees (utilisees pour le fit)
  inline int_4 NDataGood()  const {return mNDataGood;}
  //! Retourne la place maximale allouee pour les donnees
  inline int_4 NDataAlloc() const {return mNDataAlloc;}

  //! Retourne 1 si point valide, sinon 0
  inline unsigned short IsValid(int_4 i) const
                {if(i>=0 && i<mNData) return mOK[i]; else return 0;}
  //! Retourne ``true'' si il y a des erreurs sur les variables d'abscisse, ``false'' sinon.
  inline bool HasXErrors() {if(mErrXP) return true; else return false;}

  //! Retourne l'abscisse pour 1 dimension (y=f(x)) donnee I
  inline double X1(int_4 i) const
                {if(i>=0 && i<mNData) return mXP[i]; else return 0.;}
  //! Retourne la 1er abscisse (X) pour (v=f(x,y,z,...)) donnee I
  inline double X(int_4 i) const
                {if(i>=0 && i<mNData) return mXP[i*mNVar]; else return 0.;}
  //! Retourne la 2sd abscisse (Y) pour (v=f(x,y,z,...)) donnee I
  inline double Y(int_4 i) const
                {if(i>=0 && i<mNData && 1<mNVar) return mXP[i*mNVar+1]; else return 0.;}
  //! etourne la 3ieme abscisse (Z) pour (v=f(x,y,z,...)) donnee I
  inline double Z(int_4 i) const
                {if(i>=0 && i<mNData && 2<mNVar) return mXP[i*mNVar+2]; else return 0.;}
  //! Retourne la Jieme abscisse (Xj) pour (v=f(x0,x1,x2,...)) donnee I
  inline double Absc(int_4 j,int_4 i) const
                {if(i>=0 && i<mNData && j<mNVar)return mXP[i*mNVar+j]; else return 0.;}
  //! Retourne la valeur de la Ieme donnee
  inline double Val(int_4 i) const
                {if(i>=0 && i<mNData) return mF[i]; else return 0.;}

  //! Retourne l'erreur (dx) sur l'abscisse pour 1 dimension (y=f(x)) donnee I
  inline double EX1(int_4 i) const
                {if(mErrXP && i>=0 && i<mNData) return mErrXP[i]; else return 0.;}
  //! Retourne l'erreur (dx) sur la 1er abscisse (X) pour (v=f(x,y,z,...)) donnee I
  inline double EX(int_4 i) const
                {if(mErrXP && i>=0 && i<mNData) return mErrXP[i*mNVar]; else return 0.;}
  //! Retourne l'erreur (dy) sur la 2sd abscisse (Y) pour (v=f(x,y,z,...)) donnee I
  inline double EY(int_4 i) const
                {if(mErrXP && i>=0 && i<mNData && 1<mNVar) return mErrXP[i*mNVar+1];
                 else return 0.;}
  //! Retourne l'erreur (dz) sur la 3ieme abscisse (Z) pour (v=f(x,y,z,...)) donnee I
  inline double EZ(int_4 i) const
                {if(mErrXP && i>=0 && i<mNData && 2<mNVar) return mErrXP[i*mNVar+2];
                 else return 0.;}
  //! Retourne l'erreur (dxj) sur la Jieme abscisse (Xj) pour (v=f(x0,x1,x2,...)) donnee I
  inline double EAbsc(int_4 j,int_4 i) const
                {if(mErrXP && i>=0 && i<mNData && j<mNVar) return mErrXP[i*mNVar+j];
                 else return 0.;}
  //! Retourne l'erreur de la Ieme donnee
  inline double EVal(int_4 i) const
                {if(i>=0 && i<mNData) return mErr[i]; else return 0.;}

  r_8* GetVec(int_4 n, r_8* ret=NULL) const;
  r_4* GetVecR4(int_4 n, r_4* ret=NULL) const;
  int_4 GetMnMx(int_4 var,int_4& imin,int_4& imax) const;
  int_4 GetMnMx(int_4 var,double& min,double& max) const;
  int_4 GetMeanSigma(int_4 var,double& mean,double& sigma,double min=1.,double max=-1.) const;
  int_4 GetMoMeMed(int_4 var,double& mode,double& mean,double& median,
                 double min=1.,double max=-1.,double coeff=0.8) const;
  int_4 GetMode(int_4 var,double& mode,double min=1.,double max=-1.,double coeff=0.8) const;
  double PolFit(int_4 varx,Poly& pol,int_4 degre,bool ey=true,double xc=0.) const;
  double PolFit(int_4 varx,int_4 vary,Poly2& pol,int_4 degre1,int_4 degre2=-1,bool ez=true
               ,double xc=0.,double yc=0.) const;
  GeneralFitData FitResidus(GeneralFit& gfit) const;
  GeneralFitData FitFunction(GeneralFit& gfit) const;

// Declaration de l interface NTuple
  virtual sa_size_t     NbLines() const;
  virtual sa_size_t     NbColumns() const;
  virtual r_8 *         GetLineD(sa_size_t n) const;
  virtual r_8           GetCell(sa_size_t n, sa_size_t k) const;
  virtual r_8           GetCell(sa_size_t n, string const & nom) const;
  virtual void          GetMinMax(sa_size_t k, double& min, double& max) const; 
  virtual void          GetMinMax(string const & nom, double& min, double& max) const;
  virtual sa_size_t     ColumnIndex(string const & nom) const;
  virtual string        ColumnName(sa_size_t k) const;
  virtual string        VarList_C(const char* nomx=NULL) const;

  enum FitErrType {DefaultError,ConstantError,SqrtError,ProporError};
  static double ComputeError(double val,double err,FitErrType errtype=DefaultError
                            ,double errscale=1.,double errmin=0.,bool nozero=false);

protected:
  void Alloc(uint_4 nVar, uint_4 nDataAlloc, uint_2 errx=0);
  void Delete();

  int_4   mNVar;        //!< number of variables
  int_4   mNDataAlloc;  //!< space allocated for datas
  int_4   mNData;       //!< number of datas
  int_4   mNDataGood;   //!< number of good datas (not killed)
  uint_2  mOk_EXP;      //!< true if errors on variables X,Y,..
  r_8* mXP;             //!< x0 y0 z0, x1 y1 z1, ..., xn yn zn 
  r_8* mErrXP;          //!< Ex0 Ey0 Ez0, Ex1 Ey1 Ez1, ..., Exn Eyn Ezn 
  r_8* mF;              //!< F0, F1, ...., Fn
  r_8* mErr;            //!< EF0, EF1, ...., EFn
  uint_2* mOK;          //!< 1 si pt valid, 0 sinon
  r_8* BuffVar;         //!< dummy buffer
  r_4* BuffVarR4;       //!< dummy buffer
};

inline ostream& operator << (ostream& s, GeneralFitData const & g)
                         {g.Show(s); return(s);}

/////////////////////////////////////////////////////////////////////////
// Classe pour la gestion de persistance

//! For persistance management
inline POutPersist& operator << (POutPersist& os, GeneralFitData & obj)
{ ObjFileIO<GeneralFitData> fio(&obj);  fio.Write(os);  return(os); }
//! For persistance management
inline PInPersist& operator >> (PInPersist& is, GeneralFitData & obj)
{ ObjFileIO<GeneralFitData> fio(&obj); is.SkipToNextObject(); fio.Read(is); return(is); }

// Classe pour la gestion de persistance
// ObjFileIO<GeneralFitData>

} // Fin du namespace

#endif