Changeset 514 in Sophya


Ignore:
Timestamp:
Oct 25, 1999, 6:43:04 PM (26 years ago)
Author:
ansari
Message:

elimination des OVector/OMatrix au profit des TVector/TMatrix cmv 25/10/99

Location:
trunk/SophyaLib
Files:
33 edited

Legend:

Unmodified
Added
Removed

  • trunk/SophyaLib/NTools/difeq.cc

    r508 r514  
    106106//++
    107107DiffEqSolver&
    108 DiffEqSolver::StartV(OVector const& yi, double t)
     108DiffEqSolver::StartV(Vector const& yi, double t)
    109109//
    110110//      Spécifie le point de départ de l'intégration de l'équadif.
    … …  
    169169
    170170//++
    171 // virtual void DiffEqSolver::SolveArr(OMatrix&  y, double* t, double tf, int n)=0
     171// virtual void DiffEqSolver::SolveArr(Matrix&  y, double* t, double tf, int n)=0
    172172//      Lance la résolution de l'équadif, jusqu'à la valeur
    173173//      tf du temps. N valeurs intermédiaires sont retournées dans
    … …  
    189189//++
    190190void
    191 DiffEqSolver::SolveV(OVector& yf, double tf)
     191DiffEqSolver::SolveV(Vector& yf, double tf)
    192192//
    193193//      Lance la résolution de l'équadif, jusqu'à la valeur
    … …  
    201201{
    202202  double t;
    203   OMatrix m(1, mFunc->NFuncReal());
     203  Matrix m(1, mFunc->NFuncReal());
    204204  SolveArr(m, &t, tf, 1);
    205205  yf.Realloc(mFunc->NFunc());
    … …  
    219219  ASSERT(mFunc->NFunc() == 1);
    220220  double t;
    221   OMatrix m(1,mFunc->NFuncReal());
     221  Matrix m(1,mFunc->NFuncReal());
    222222  SolveArr(m, &t, tf, 1);
    223223  yf = m(0,0);
    … …  
    238238{
    239239  double t;
    240   OMatrix m(1, mFunc->NFuncReal());
     240  Matrix m(1, mFunc->NFuncReal());
    241241  SolveArr(m, &t, tf, 1);
    242242  for (int i=0; i<mFunc->NFunc(); i++)
    … …  
    258258{
    259259  ASSERT(mFunc->NFunc() == 1);
    260   OMatrix m(n, mFunc->NFuncReal());
     260  Matrix m(n, mFunc->NFuncReal());
    261261  SolveArr(m, t, tf, n);
    262262  for (int i=0; i<n; i++)
    … …  
    281281//--
    282282{
    283    OMatrix m(n, mFunc->NFuncReal());
     283   Matrix m(n, mFunc->NFuncReal());
    284284   SolveArr(m, t, tf, n);
    285285   for (int i=0; i<n; i++)
    … …  
    318318
    319319//++
    320 //  virtual void ComputeV(OVector& fpi, OVector const& fi)
     320//  virtual void ComputeV(Vector& fpi, Vector const& fi)
    321321//      Calcule les valeurs des dérivées fpi à partir des valeurs
    322322//      des fonctions fi. A redéfinir.
    … …  
    340340
    341341//++
    342 // virtual void AdjustStart(OVector& start, double tstart)
     342// virtual void AdjustStart(Vector& start, double tstart)
    343343//      Pour ajuster le vecteur de départ quand il y a des
    344344//      fonctions à usage interne...
    … …  
    405405
    406406void
    407 DiffEqFcnT1::ComputeV(OVector& fpi, OVector const& fi)
     407DiffEqFcnT1::ComputeV(Vector& fpi, Vector const& fi)
    408408{
    409409  fpi(0) = (*mFcn)(fi(0), fi(1));
    … …  
    412412
    413413void
    414 DiffEqFcnT1::AdjustStart(OVector& start, double tstart)
     414DiffEqFcnT1::AdjustStart(Vector& start, double tstart)
    415415{
    416416  start.Realloc(2);
    … …  
    446446
    447447void
    448 DiffEqFcn2::ComputeV(OVector& fpi, OVector const& fi)
     448DiffEqFcn2::ComputeV(Vector& fpi, Vector const& fi)
    449449{
    450450  fpi(0) = fi(1);
    … …  
    481481
    482482void
    483 DiffEqFcnT2::ComputeV(OVector& fpi, OVector const& fi)
     483DiffEqFcnT2::ComputeV(Vector& fpi, Vector const& fi)
    484484{
    485485  fpi(0) = fi(1);
    … …  
    489489
    490490void
    491 DiffEqFcnT2::AdjustStart(OVector& start, double tstart)
     491DiffEqFcnT2::AdjustStart(Vector& start, double tstart)
    492492{
    493493  start.Realloc(3);
    … …  
    505505//      un vecteur de dimension 3.
    506506//      On fournit une fonction de type
    507 //|     typedef void(*DIFEQFCNV)(OVector& y2, OVector const& y1,
    508 //|                              OVector const& y);
     507//|     typedef void(*DIFEQFCNV)(Vector& y2, Vector const& y1,
     508//|                              Vector const& y);
    509509//      qui retourne y'' en fonction de y' et y.
    510510//      Note : le système résolu est alors en fait
    … …  
    527527
    528528void
    529 DiffEqFcnV::ComputeV(OVector& fpi, OVector const& fi)
     529DiffEqFcnV::ComputeV(Vector& fpi, Vector const& fi)
    530530{
    531531  fpi(0) = fi(3);
    … …  
    582582
    583583void
    584 RK4DiffEq::SolveArr(OMatrix& y, double* t, double tf, int n)
     584RK4DiffEq::SolveArr(Matrix& y, double* t, double tf, int n)
    585585{
    586586  //TIMEF;
    … …  
    592592   double dx = (tf - mXStart)/(n*nStep);
    593593   
    594    OVector yt = mYStart;
     594   Vector yt = mYStart;
    595595
    596596   k1.Realloc(mFunc->NFuncReal());
    … …  
    610610
    611611void
    612 RK4DiffEq::RKStep(OVector& newY, OVector const& y0, double dt)
     612RK4DiffEq::RKStep(Vector& newY, Vector const& y0, double dt)
    613613{
    614614    mFunc->ComputeV(k1, y0);
    615615    k1 *= dt;
    616     mFunc->ComputeV(k2, y0 + k1/2);
     616    mFunc->ComputeV(k2, y0 + k1/2.);
    617617    k2 *= dt;
    618     mFunc->ComputeV(k3, y0 + k2/2);
     618    mFunc->ComputeV(k3, y0 + k2/2.);
    619619    k3 *= dt;
    620620    mFunc->ComputeV(k4, y0 + k3);
    621621    k4 *= dt;
    622622   
    623     newY = y0 + (k1 + k2*2. + k3*2. + k4)/6;
    624 }
    625 
     623    newY = y0 + (k1 + k2*2. + k3*2. + k4)/6.;
     624}
     625

  • trunk/SophyaLib/NTools/difeq.h

    r508 r514  
    55#include "machdefs.h"
    66#include "pexceptions.h"
    7 #include "cvector.h"
    8 
    9 namespace PlanckDPC {class GeneralFunction;}
     7#include "tvector.h"
     8
     9namespace PlanckDPC {
     10
     11class GeneralFunction;
    1012
    1113// <summary> fonction pour equadifs </summary>
    … …  
    2830 
    2931  // Calcule les valeurs des derivees fpi a partir des valeurs des fonctions fi
    30   virtual void ComputeV(OVector& fpi, OVector const& fi)
     32  virtual void ComputeV(Vector& fpi, Vector const& fi)
    3133       { Compute(fpi(0), fi(0)); }
    3234
    … …  
    4446  // Pour ajuster vecteur de depart quand il y a des fonctions a usage
    4547  // interne...
    46   virtual void AdjustStart(OVector& /*start*/, double /*tstart*/)
     48  virtual void AdjustStart(Vector& /*start*/, double /*tstart*/)
    4749    {}   
    4850protected:
    … …  
    8789  // Implementation de Compute qui va utiliser la fonction fournie
    8890  // au constructeur.
    89   virtual void ComputeV(OVector& fpi, OVector const& fi);
     91  virtual void ComputeV(Vector& fpi, Vector const& fi);
    9092
    9193  // Implementation de AdjustStart qui gere la fonction a usage interne.
    92   virtual void AdjustStart(OVector& start, double tstart);
     94  virtual void AdjustStart(Vector& start, double tstart);
    9395protected:
    9496  DIFEQFCNT1 mFcn;
    … …  
    102104  DiffEqFcn2(DIFEQFCN2);
    103105
    104   virtual void ComputeV(OVector& fpi, OVector const& fi);
     106  virtual void ComputeV(Vector& fpi, Vector const& fi);
    105107protected:
    106108  DIFEQFCN2 mFcn;
    … …  
    121123  // Implementation de Compute qui va utiliser la fonction fournie
    122124  // au constructeur.
    123   virtual void ComputeV(OVector& fpi, OVector const& fi);
     125  virtual void ComputeV(Vector& fpi, Vector const& fi);
    124126
    125127  // Implementation de AdjustStart qui gere la fonction a usage interne.
    126   virtual void AdjustStart(OVector& start, double tstart);
     128  virtual void AdjustStart(Vector& start, double tstart);
    127129protected:
    128130  DIFEQFCNT2 mFcn;
    … …  
    130132
    131133// Cas y'' = f(y',y) avec des 3-vecteurs
    132 typedef void(*DIFEQFCNV)(OVector&, OVector const&, OVector const&);
     134typedef void(*DIFEQFCNV)(Vector&, Vector const&, Vector const&);
    133135
    134136// <summary> y'' = f(y',y,t) </summary>
    135137// Cas y'' = f(y',y,t), on fournit la fonction f, sous la forme
    136 // double f(OVector), et ca construit la bonne DiffEqFunction
     138// double f(Vector), et ca construit la bonne DiffEqFunction
    137139class DiffEqFcnV : public DiffEqFunction {
    138140public:
    139   // Constructeur, on fournit une fonction (OVector)->double
     141  // Constructeur, on fournit une fonction (Vector)->double
    140142  // qui donne y'' en fonction du vecteur (t, y, y')
    141143  DiffEqFcnV(DIFEQFCNV);
    … …  
    143145  // Implementation de Compute qui va utiliser la fonction fournie
    144146  // au constructeur.
    145   virtual void ComputeV(OVector& fpi, OVector const& fi);
     147  virtual void ComputeV(Vector& fpi, Vector const& fi);
    146148protected:
    147149  DIFEQFCNV mFcn;
    148   OVector tmp1, tmp2, tmp3;
     150  Vector tmp1, tmp2, tmp3;
    149151};
    150152
    … …  
    176178  // Change les conditions initiales. Notation chainee possible.
    177179  // <group>
    178   DiffEqSolver& StartV(OVector const& yi, double t);
     180  DiffEqSolver& StartV(Vector const& yi, double t);
    179181  // si NFunc == 1
    180182  DiffEqSolver& Start1(double        yi, double t);
    … …  
    184186  // Lance la resolution, avec ou sans conservation de n valeurs intermediaires
    185187  // <group>
    186   virtual void SolveV(OVector& yf, double tf);
     188  virtual void SolveV(Vector& yf, double tf);
    187189  // si NFunc == 1
    188190  virtual void Solve1(double& yf, double tf); 
    189191  virtual void Solve(double* yf, double tf);
    190   virtual void SolveArr(OMatrix&  y, double* t, double tf, int n)=0;
     192  virtual void SolveArr(Matrix&  y, double* t, double tf, int n)=0;
    191193  // si NFunc == 1
    192194  virtual void SolveArr1(double*  y, double* t, double tf, int n);   
    … …  
    198200  bool mOwnFunc;
    199201
    200   OVector  mYStart;
     202  Vector  mYStart;
    201203  double  mXStart;
    202204  double  mStep;
    … …  
    215217
    216218  // Implementation de RK4
    217   virtual void SolveArr(OMatrix& y, double* t, double tf, int n);
     219  virtual void SolveArr(Matrix& y, double* t, double tf, int n);
    218220
    219221protected:
    220222  // Un pas RK4
    221   void RKStep(OVector& newY, OVector const& y0, double dt);
     223  void RKStep(Vector& newY, Vector const& y0, double dt);
    222224  // Vecteurs utilises en interne, pour ne pas les reallouer.
    223   OVector k1, k2, k3, k4;
    224 };
    225 
    226 
     225  Vector k1, k2, k3, k4;
     226};
     227
     228
     229} // Fin du namespace
    227230
    228231#endif

  • trunk/SophyaLib/NTools/dvlist.cc

    r404 r514  
    3737// PPersist
    3838//--
    39 
    40 
    41 using namespace PlanckDPC;
    4239
    4340char MuTyV::myStrBuf[64];   // Declare static ds le .h

  • trunk/SophyaLib/NTools/fct1dfit.cc

    r307 r514  
    88#include "nbconst.h"
    99#include "tabmath.h"
    10 
    11 using namespace PlanckDPC;
    1210
    1311//define EXPO exp

  • trunk/SophyaLib/NTools/fct2dfit.cc

    r490 r514  
    1616#define EXPO tabFExp
    1717#define MINEXPM (100.)
    18 
    19 using namespace PlanckDPC;
    2018
    2119//================================================================

  • trunk/SophyaLib/NTools/fftserver.cc

    r508 r514  
    156156}
    157157
    158 void FFTServer::fftf(OVector& in, OVector& out)
     158void FFTServer::fftf(Vector& in, Vector& out)
    159159{
    160160  int l = in.NElts();
    161 /* ----- Si c'etait un OVector<float> on aurait ecrit comme ca
    162    Pour le moment OVector est double, on passe donc par un tableau
     161/* ----- Si c'etait un Vector<float> on aurait ecrit comme ca
     162   Pour le moment Vector est double, on passe donc par un tableau
    163163   intermediare
    164164// La transformee sur le tableau de flaot se fait en place,
    … …  
    175175}
    176176
    177 void FFTServer::fftb(OVector& in, OVector& out)
     177void FFTServer::fftb(Vector& in, Vector& out)
    178178{
    179179  int l = in.NElts();
    180 /* ----- Si c'etait un OVector<float> on aurait ecrit comme ca
    181    Pour le moment OVector est double, on passe donc par un tableau
     180/* ----- Si c'etait un Vector<float> on aurait ecrit comme ca
     181   Pour le moment Vector est double, on passe donc par un tableau
    182182   intermediare
    183183// La transformee sur le tableau de flaot se fait en place,

  • trunk/SophyaLib/NTools/fftserver.h

    r508 r514  
    33
    44#include <complex>
    5 #include "cvector.h"
     5#include "tvector.h"
    66
    77class FFTServer{
    … …  
    1717  virtual void fftf(int l, complex<double>* inout);
    1818  virtual void fftb(int l, complex<double>* inout);
    19   virtual void fftf(OVector& in, OVector& out);
    20   virtual void fftb(OVector& in, OVector& out);
     19  virtual void fftf(Vector& in, Vector& out);
     20  virtual void fftb(Vector& in, Vector& out);
    2121
    2222 protected:
    … …  
    9494  \param inout input array /output backward FFT(original array destroyed)
    9595*/
    96 /*!\fn  virtual void FFTServer::fftf(OVector& in, OVector& out)
     96/*!\fn  virtual void FFTServer::fftf(Vector& in, Vector& out)
    9797  \param in input array
    9898  \param out forward FFT
    9999*/
    100 /*! \fn virtual void FFTServer::fftb(OVector& in, OVector& out)
     100/*! \fn virtual void FFTServer::fftb(Vector& in, Vector& out)
    101101  \param in input array
    102102  \param out backward FFT

  • trunk/SophyaLib/NTools/generaldata.cc

    r508 r514  
    1414#include "pexceptions.h"
    1515#include "objfio.h"
    16 
    17 using namespace PlanckDPC;
    1816
    1917//================================================================
    … …  
    748746if(varx<0 || varx>=mNVar) return -2.;
    749747if(mNDataGood<=0) return -4.;
    750 OVector x(mNDataGood);
    751 OVector y(mNDataGood);
    752 OVector ey2(1);
     748Vector x(mNDataGood);
     749Vector y(mNDataGood);
     750Vector ey2(1);
    753751if(ey) ey2.Realloc(mNDataGood,true);
    754752int ntest = 0;
    … …  
    763761double res = 0.;
    764762if(ey) {
    765   OVector errcoef(1);
     763  Vector errcoef(1);
    766764  res = pol.Fit(x,y,ey2,degre,errcoef);
    767765} else {
    … …  
    797795if(vary<0 || vary>=mNVar || vary==varx) return -3.;
    798796if(mNDataGood<=0) return -4.;
    799 OVector x(mNDataGood);
    800 OVector y(mNDataGood);
    801 OVector z(mNDataGood);
    802 OVector ez2(1);
     797Vector x(mNDataGood);
     798Vector y(mNDataGood);
     799Vector z(mNDataGood);
     800Vector ez2(1);
    803801if(ez) ez2.Realloc(mNDataGood,true);
    804802int ntest = 0;
    … …  
    814812double res = 0.;
    815813if(ez) {
    816   OVector errcoef(1);
     814  Vector errcoef(1);
    817815  if(degre2>0) res = pol.Fit(x,y,z,ez2,degre1,degre2,errcoef);
    818816  else         res = pol.Fit(x,y,z,ez2,degre1,errcoef);

  • trunk/SophyaLib/NTools/generalfit.cc

    r508 r514  
    99#include "pexceptions.h"
    1010#include "generalfit.h"
    11 #include "cvector.h"
    1211
    1312#define EPS_FIT_MIN 1.e-8
    14 
    15 using namespace PlanckDPC;
    1613
    1714//================================================================
    … …  
    874871
    875872//++
    876 OVector GeneralFit::GetParm()
     873Vector GeneralFit::GetParm()
    877874//
    878875//      Retourne les valeurs des parametres dans un vecteur.
    … …  
    10961093{
    10971094 volatile double oldChi2;
    1098  OMatrix COVAR(mNPar,mNPar);
    1099  OVector DA(mNPar);
    1100  OVector dparam(mNPar);
    1101  OVector paramTry(mNPar);
    1102  OVector param_tr(mNPar);
    1103  OVector paramTry_tr(mNPar);
    1104  OVector step_tr(mNPar);
     1095 Matrix COVAR(mNPar,mNPar);
     1096 Vector DA(mNPar);
     1097 Vector dparam(mNPar);
     1098 Vector paramTry(mNPar);
     1099 Vector param_tr(mNPar);
     1100 Vector paramTry_tr(mNPar);
     1101 Vector step_tr(mNPar);
    11051102 nStop = nStopL = nStep = 0;
    11061103 Chi2 = oldChi2 = 0.;
    … …  
    15151512
    15161513//////////////////////////////////////////////////////////////////////
    1517 void GeneralFit::write_in_step(double ci2,OVector& par)
     1514void GeneralFit::write_in_step(double ci2,Vector& par)
    15181515{
    15191516if(FileStep==NULL) return;
    … …  
    15241521
    15251522//////////////////////////////////////////////////////////////////////
    1526 void GeneralFit::TryFunc(OVector& par,OVector& par_tr)
     1523void GeneralFit::TryFunc(Vector& par,Vector& par_tr)
    15271524{
    15281525  BETA_Try = 0;
    15291526  ATGA_Try = 0;
    15301527  Chi2  = 0;
    1531   OVector deriv(mNPar);
    1532   OVector derivtr(mNPar);
     1528  Vector deriv(mNPar);
     1529  Vector derivtr(mNPar);
    15331530  double result;
    15341531
    … …  
    15641561
    15651562//////////////////////////////////////////////////////////////////////
    1566 void GeneralFit::TryXi2(OVector& par,OVector& par_tr)
     1563void GeneralFit::TryXi2(Vector& par,Vector& par_tr)
    15671564{
    15681565  double c, *parloc;
    … …  
    16641661
    16651662//////////////////////////////////////////////////////////////////////
    1666 OVector GeneralFit::p_vers_tr(OVector const& p)
    1667 {
    1668  OVector tr(p);
     1663Vector GeneralFit::p_vers_tr(Vector const& p)
     1664{
     1665 Vector tr(p);
    16691666 for(int i=0;i<mNPar;i++) {
    16701667   if( fixParam[i] || ! boundParam[i] ) continue;
    … …  
    16751672
    16761673//////////////////////////////////////////////////////////////////////
    1677 void GeneralFit::p_vers_tr(OVector const& p,OVector& tr)
     1674void GeneralFit::p_vers_tr(Vector const& p,Vector& tr)
    16781675{
    16791676 for(int i=0;i<mNPar;i++) {
    … …  
    16951692
    16961693//////////////////////////////////////////////////////////////////////
    1697 OVector GeneralFit::tr_vers_p(OVector const& tr)
    1698 {
    1699  OVector p(tr);
     1694Vector GeneralFit::tr_vers_p(Vector const& tr)
     1695{
     1696 Vector p(tr);
    17001697 for(int i=0;i<mNPar;i++) {
    17011698   if( fixParam[i] || ! boundParam[i] ) continue;
    … …  
    17061703
    17071704//////////////////////////////////////////////////////////////////////
    1708 void GeneralFit::tr_vers_p(OVector const& tr,OVector& p)
     1705void GeneralFit::tr_vers_p(Vector const& tr,Vector& p)
    17091706{
    17101707 for(int i=0;i<mNPar;i++) {
    … …  
    17271724
    17281725//////////////////////////////////////////////////////////////////////
    1729 OVector GeneralFit::dp_vers_dtr(OVector const& dp,OVector const& tr)
    1730 {
    1731  OVector dtr(dp);
     1726Vector GeneralFit::dp_vers_dtr(Vector const& dp,Vector const& tr)
     1727{
     1728 Vector dtr(dp);
    17321729 for(int i=0;i<mNPar;i++) {
    17331730   if( fixParam[i] || ! boundParam[i] ) continue;
    … …  
    17381735
    17391736//////////////////////////////////////////////////////////////////////
    1740 void GeneralFit::dp_vers_dtr(OVector const& dp,OVector const& tr,OVector& dtr)
     1737void GeneralFit::dp_vers_dtr(Vector const& dp,Vector const& tr,Vector& dtr)
    17411738{
    17421739 for(int i=0;i<mNPar;i++) {
    … …  
    17591756
    17601757//////////////////////////////////////////////////////////////////////
    1761 OVector GeneralFit::dtr_vers_dp(OVector const& dtr,OVector const& tr)
    1762 {
    1763  OVector dp(dtr);
     1758Vector GeneralFit::dtr_vers_dp(Vector const& dtr,Vector const& tr)
     1759{
     1760 Vector dp(dtr);
    17641761 for(int i=0;i<mNPar;i++) {
    17651762   if( fixParam[i] || ! boundParam[i] ) continue;
    … …  
    17711768//////////////////////////////////////////////////////////////////////
    17721769// inline fonction pour aller + vite dans le try()
    1773 //void GeneralFit::dtr_vers_dp(OVector const& dtr,OVector const& tr,OVector& dp)
    1774 
    1775 //////////////////////////////////////////////////////////////////////
    1776 int GeneralFit::put_in_limits_for_deriv(OVector const& p,OVector& dp,double dist)
     1770//void GeneralFit::dtr_vers_dp(Vector const& dtr,Vector const& tr,Vector& dp)
     1771
     1772//////////////////////////////////////////////////////////////////////
     1773int GeneralFit::put_in_limits_for_deriv(Vector const& p,Vector& dp,double dist)
    17771774// 1-/ Redefinit dp pour qu'il soit superieur a minStepDeriv
    17781775// 2-/ Redefinit dp pour que p+/-dp reste dans les limites (parametre borne)

  • trunk/SophyaLib/NTools/generalfit.h

    r508 r514  
    44
    55#include "pexceptions.h"
    6 #include "matrix.h"
    7 #include "cvector.h"
     6#include "tvector.h"
    87#include "generaldata.h"
    98
    … …  
    124123
    125124  double          GetParm(int);
    126   OVector          GetParm();
     125  Vector          GetParm();
    127126  double          GetParmErr(int);
    128127  double          GetCoVar(int,int);
    … …  
    168167  GeneralFitData*   mData;
    169168
    170   OVector          Param;
    171   OVector          errParam;
    172   OVector          stepParam;
    173   OVector          minParam;
    174   OVector          maxParam;
    175   OVector          minStepDeriv;
    176   OVector          Eps;
     169  Vector          Param;
     170  Vector          errParam;
     171  Vector          stepParam;
     172  Vector          minParam;
     173  Vector          maxParam;
     174  Vector          minStepDeriv;
     175  Vector          Eps;
    177176  unsigned short int* fixParam;
    178177  unsigned short int* boundParam;
    … …  
    188187  FILE            *FileStep;
    189188 
    190   OMatrix          ATGA;
    191   OVector          BETA;
    192   OMatrix          ATGA_Try;
    193   OVector          BETA_Try;
    194   OVector          C;
    195   OVector          D;
     189  Matrix          ATGA;
     190  Vector          BETA;
     191  Matrix          ATGA_Try;
     192  Vector          BETA_Try;
     193  Vector          C;
     194  Vector          D;
    196195 
    197196  double          Chi2;
    … …  
    202201 
    203202  // Fonctions privees
    204   void       write_in_step(double ci2,OVector& par);
     203  void       write_in_step(double ci2,Vector& par);
    205204  void       General_Init(void);
    206   void       TryFunc(OVector& par,OVector& par_tr);
    207   void       TryXi2(OVector& par,OVector& par_tr);
     205  void       TryFunc(Vector& par,Vector& par_tr);
     206  void       TryXi2(Vector& par,Vector& par_tr);
    208207  void       CheckSanity();
    209208  void       Set_Bound_C_D(int i);
    210209  void       Set_Bound_C_D();
    211210  double     p_vers_tr(int i,double p);
    212   OVector     p_vers_tr(OVector const& p);
    213   void       p_vers_tr(OVector const& p,OVector& tr);
     211  Vector     p_vers_tr(Vector const& p);
     212  void       p_vers_tr(Vector const& p,Vector& tr);
    214213  double     tr_vers_p(int i,double tr);
    215   OVector     tr_vers_p(OVector const& tr);
    216   void       tr_vers_p(OVector const& tr,OVector& p);
     214  Vector     tr_vers_p(Vector const& tr);
     215  void       tr_vers_p(Vector const& tr,Vector& p);
    217216  double     c_dp_vers_dtr(int i,double tr);
    218   OVector       dp_vers_dtr(OVector const& dp,OVector const& tr);
    219   void         dp_vers_dtr(OVector const& dp,OVector const& tr,OVector& dtr);
     217  Vector       dp_vers_dtr(Vector const& dp,Vector const& tr);
     218  void         dp_vers_dtr(Vector const& dp,Vector const& tr,Vector& dtr);
    220219  double     c_dtr_vers_dp(int i,double tr);
    221   OVector       dtr_vers_dp(OVector const& dtr,OVector const& tr);
    222   int        put_in_limits_for_deriv(OVector const& p,OVector& dp,double dist=0.66);
    223   inline void dtr_vers_dp(OVector const& dtr,OVector const& tr,OVector& dp)
     220  Vector       dtr_vers_dp(Vector const& dtr,Vector const& tr);
     221  int        put_in_limits_for_deriv(Vector const& p,Vector& dp,double dist=0.66);
     222  inline void dtr_vers_dp(Vector const& dtr,Vector const& tr,Vector& dp)
    224223              {  for(int i=0;i<mNPar;i++)
    225224                   { if( fixParam[i] ) continue;

  • trunk/SophyaLib/NTools/hisprof.cc

    r508 r514  
    55#include "hisprof.h"
    66#include "perrors.h"
    7 #include "cvector.h"
    8 
    9 using namespace PlanckDPC;
    107
    118//++
    … …  
    260257//--
    261258//++
    262 // inline void GetMean(OVector& v)
     259// inline void GetMean(Vector& v)
    263260//      Retourne le contenu de la moyenne dans le vecteur v
    264261//--
    265262//++
    266 // inline void GetError2(OVector& v)
     263// inline void GetError2(Vector& v)
    267264//      Retourne le contenu au carre de la dispersion/erreur dans le vecteur v
    268265//--

  • trunk/SophyaLib/NTools/hisprof.h

    r508 r514  
    44#include <stdio.h>
    55#include "peida.h"
    6 #include "cvector.h"
     6#include "tvector.h"
    77#include "ppersist.h"
    88#include "histos.h"
    … …  
    2929  // Acces a l information
    3030  inline Histo GetHisto() {if(!Ok) UpdateHisto(); return (Histo) *this;}
    31   inline void GetMean(OVector& v) {if(!Ok) UpdateHisto(); Histo::GetValue(v);}
    32   inline void GetError2(OVector& v) {if(!Ok) UpdateHisto(); Histo::GetError2(v);}
     31  inline void GetMean(Vector& v) {if(!Ok) UpdateHisto(); Histo::GetValue(v);}
     32  inline void GetError2(Vector& v) {if(!Ok) UpdateHisto(); Histo::GetError2(v);}
    3333  inline float operator()(int i) const {if(!Ok) UpdateHisto(); return data[i];}
    3434  inline float Error2(int i) const {if(!Ok) UpdateHisto(); return (float) err2[i];}

  • trunk/SophyaLib/NTools/histos.cc

    r508 r514  
    11//
    2 // $Id: histos.cc,v 1.5 1999-10-25 10:36:07 ansari Exp $
     2// $Id: histos.cc,v 1.6 1999-10-25 16:39:57 ansari Exp $
    33//
    44
    … …  
    99#include "histos.h"
    1010#include "perrors.h"
    11 #include "cvector.h"
    1211#include "poly.h"
    1312#include "strutil.h"
    1413#include "generalfit.h"
    15 
    16 using namespace PlanckDPC;
    1714
    1815//++
    … …  
    404401/********* Methode *********/
    405402//++
    406 void Histo::GetAbsc(OVector &v)
     403void Histo::GetAbsc(Vector &v)
    407404//
    408405//      Remplissage d'un tableau avec la valeur des abscisses
    … …  
    415412
    416413//++
    417 void Histo::GetValue(OVector &v)
     414void Histo::GetValue(Vector &v)
    418415//
    419416//      Remplissage d'un tableau avec la valeur du contenu
    … …  
    426423
    427424//++
    428 void Histo::GetError2(OVector &v)
     425void Histo::GetError2(Vector &v)
    429426//
    430427//      Remplissage d'un tableau avec la valeur des erreurs au carre
    … …  
    438435
    439436//++
    440 void Histo::GetError(OVector &v)
     437void Histo::GetError(Vector &v)
    441438//
    442439//      Remplissage d'un tableau avec la valeur des erreurs
    … …  
    451448/********* Methode *********/
    452449//++
    453 void Histo::PutValue(OVector &v, int ierr)
     450void Histo::PutValue(Vector &v, int ierr)
    454451//
    455452//      Remplissage du contenu de l'histo avec les valeurs d'un vecteur
    … …  
    465462
    466463//++
    467 void Histo::PutValueAdd(OVector &v, int ierr)
     464void Histo::PutValueAdd(Vector &v, int ierr)
    468465//
    469466//      Addition du contenu de l'histo avec les valeurs d'un vecteur
    … …  
    479476
    480477//++
    481 void Histo::PutError2(OVector &v)
     478void Histo::PutError2(Vector &v)
    482479//
    483480//      Remplissage des erreurs au carre de l'histo avec les valeurs d'un vecteur
    … …  
    491488
    492489//++
    493 void Histo::PutError2Add(OVector &v)
     490void Histo::PutError2Add(Vector &v)
    494491//
    495492//      Addition des erreurs au carre de l'histo avec les valeurs d'un vecteur
    … …  
    503500
    504501//++
    505 void Histo::PutError(OVector &v)
     502void Histo::PutError(Vector &v)
    506503//
    507504//      Remplissage des erreurs de l'histo avec les valeurs d'un vecteur
    … …  
    11211118  }
    11221119
    1123   OVector xFit(nLowHigh);
    1124   OVector yFit(nLowHigh);
    1125   OVector e2Fit(nLowHigh);
    1126   OVector errcoef(1);
     1120  Vector xFit(nLowHigh);
     1121  Vector yFit(nLowHigh);
     1122  Vector e2Fit(nLowHigh);
     1123  Vector errcoef(1);
    11271124  int ii = 0;
    11281125  for (int j=iLow; j<=iHigh; j++) {
    … …  
    14441441if(f==NULL)
    14451442  throw(NullPtrError("Histo::FitResidus: NULL pointer\n"));
    1446 OVector par = gfit.GetParm();
     1443Vector par = gfit.GetParm();
    14471444Histo h(*this);
    14481445for(int i=0;i<NBins();i++) {
    … …  
    14641461if(f==NULL)
    14651462  throw(NullPtrError("Histo::FitFunction: NULL pointer\n"));
    1466 OVector par = gfit.GetParm();
     1463Vector par = gfit.GetParm();
    14671464Histo h(*this);
    14681465for(int i=0;i<NBins();i++) {

  • trunk/SophyaLib/NTools/histos.h

    r508 r514  
    11// This may look like C code, but it is really -*- C++ -*-
    22//
    3 // $Id: histos.h,v 1.4 1999-10-25 10:36:07 ansari Exp $
     3// $Id: histos.h,v 1.5 1999-10-25 16:39:57 ansari Exp $
    44//
    55
    … …  
    1111#include <stdio.h>
    1212#include "peida.h"
    13 #include "cvector.h"
     13#include "tvector.h"
    1414#include "ppersist.h"
    1515#include "anydataobj.h"
    … …  
    6666
    6767  // get/put dans/depuis un vector
    68   void GetAbsc(OVector& v);
    69   void GetValue(OVector& v);
    70   void GetError2(OVector& v);
    71   void GetError(OVector& v);
    72   void PutValue(OVector& v, int ierr=0);
    73   void PutValueAdd(OVector &v, int ierr=0);
    74   void PutError2(OVector& v);
    75   void PutError2Add(OVector& v);
    76   void PutError(OVector& v);
     68  void GetAbsc(Vector& v);
     69  void GetValue(Vector& v);
     70  void GetError2(Vector& v);
     71  void GetError(Vector& v);
     72  void PutValue(Vector& v, int ierr=0);
     73  void PutValueAdd(Vector &v, int ierr=0);
     74  void PutError2(Vector& v);
     75  void PutError2Add(Vector& v);
     76  void PutError(Vector& v);
    7777
    7878  // INLINES

  • trunk/SophyaLib/NTools/histos2.cc

    r508 r514  
    1515#include "histos2.h"
    1616#include "generalfit.h"
    17 
    18 using namespace PlanckDPC;
    1917
    2018//++
    … …  
    603601///////////////////////////////////////////////////////////////////
    604602//++
    605 void Histo2D::GetXCoor(OVector &v)
     603void Histo2D::GetXCoor(Vector &v)
    606604//
    607605//      Remplissage d'un tableau avec les valeurs des abscisses.
    … …  
    615613
    616614//++
    617 void Histo2D::GetYCoor(OVector &v)
     615void Histo2D::GetYCoor(Vector &v)
    618616//
    619617//      Remplissage d'un tableau avec les valeurs des ordonnees.
    … …  
    644642
    645643//++
    646 void Histo2D::GetValue(OMatrix &v)
     644void Histo2D::GetValue(Matrix &v)
    647645//
    648646//      Remplissage d'un tableau avec les valeurs du contenu.
    … …  
    656654
    657655//++
    658 void Histo2D::GetError2(OMatrix &v)
     656void Histo2D::GetError2(Matrix &v)
    659657//
    660658//      Remplissage d'un tableau avec les valeurs du carre des erreurs.
    … …  
    670668
    671669//++
    672 void Histo2D::GetError(OMatrix &v)
     670void Histo2D::GetError(Matrix &v)
    673671//
    674672//      Remplissage d'un tableau avec les valeurs des erreurs.
    … …  
    685683///////////////////////////////////////////////////////////////////
    686684//++
    687 void Histo2D::PutValue(OMatrix &v, int ierr)
     685void Histo2D::PutValue(Matrix &v, int ierr)
    688686//
    689687//      Remplissage du contenu de l'histo avec les valeurs d'un tableau.
    … …  
    700698
    701699//++
    702 void Histo2D::PutValueAdd(OMatrix &v, int ierr)
     700void Histo2D::PutValueAdd(Matrix &v, int ierr)
    703701//
    704702//      Addition du contenu de l'histo avec les valeurs d'un tableau.
    … …  
    715713
    716714//++
    717 void Histo2D::PutError2(OMatrix &v)
     715void Histo2D::PutError2(Matrix &v)
    718716//
    719717//      Remplissage des erreurs au carre de l'histo
    … …  
    729727
    730728//++
    731 void Histo2D::PutError2Add(OMatrix &v)
     729void Histo2D::PutError2Add(Matrix &v)
    732730//
    733731//      Addition des erreurs au carre de l'histo
    … …  
    744742
    745743//++
    746 void Histo2D::PutError(OMatrix &v)
     744void Histo2D::PutError(Matrix &v)
    747745//
    748746//      Remplissage des erreurs de l'histo avec les valeurs d'un tableau.
    … …  
    11051103if(f==NULL)
    11061104  throw(NullPtrError("Histo2D::FitResidus: NULL pointer\n"));
    1107 OVector par = gfit.GetParm();
     1105Vector par = gfit.GetParm();
    11081106Histo2D h2(*this);
    11091107for(int i=0;i<NBinX();i++) for(int j=0;j<NBinY();j++) {
    … …  
    11271125if(f==NULL)
    11281126  throw(NullPtrError("Histo2D::FitFunction: NULL pointer\n"));
    1129 OVector par = gfit.GetParm();
     1127Vector par = gfit.GetParm();
    11301128Histo2D h2(*this);
    11311129for(int i=0;i<NBinX();i++) for(int j=0;j<NBinY();j++) {

  • trunk/SophyaLib/NTools/histos2.h

    r508 r514  
    6363
    6464  // get/put dans/depuis une matrice / vector
    65   void GetXCoor(OVector& v);
    66   void GetValue(OMatrix &v);
    67   void GetYCoor(OVector& v);
    68   void GetError2(OMatrix& v);
    69   void GetError(OMatrix& v);
    70   void PutValue(OMatrix& v, int ierr=0);
    71   void PutValueAdd(OMatrix& v, int ierr=0);
    72   void PutError2(OMatrix& v);
    73   void PutError2Add(OMatrix& v);
    74   void PutError(OMatrix& v);
     65  void GetXCoor(Vector& v);
     66  void GetValue(Matrix &v);
     67  void GetYCoor(Vector& v);
     68  void GetError2(Matrix& v);
     69  void GetError(Matrix& v);
     70  void PutValue(Matrix& v, int ierr=0);
     71  void PutValueAdd(Matrix& v, int ierr=0);
     72  void PutError2(Matrix& v);
     73  void PutError2Add(Matrix& v);
     74  void PutError(Matrix& v);
    7575
    7676  // INLINES

  • trunk/SophyaLib/NTools/integ.h

    r490 r514  
    33// integ.h
    44//
    5 // $Id: integ.h,v 1.4 1999-10-21 15:25:49 ansari Exp $
     5// $Id: integ.h,v 1.5 1999-10-25 16:39:59 ansari Exp $
    66//
    77
    … …  
    1313#include <set>
    1414
    15 namespace PlanckDPC {class GeneralFunction;}
     15namespace PlanckDPC {
     16
     17class GeneralFunction;
    1618
    1719class Integrator EXC_AWARE {
    … …  
    9597};
    9698
    97 
     99} // Fin du namespace
    98100
    99101#endif

  • trunk/SophyaLib/NTools/linfit.cc

    r508 r514  
    22#include "peida.h"
    33#include "linfit.h"
    4 #include "matrix.h"
    5 #include "cvector.h"
    64
    7 double LinFit(const OVector& x, const OVector& y, int nf, double (*f)(int, double),
    8            OVector& c)
     5double LinFit(const Vector& x, const Vector& y, int nf, double (*f)(int, double),
     6           Vector& c)
    97{
    108  int n = x.NElts();
    119  if (n != y.NElts()) THROW(sizeMismatchErr);
    1210 
    13   OMatrix fx(nf, n);
     11  Matrix fx(nf, n);
    1412  for (int i=0; i<nf; i++)
    1513    for (int j=0; j<n; j++)
    … …  
    2119
    2220
    23 double LinFit(const OMatrix& fx, const OVector& y, OVector& c)
     21double LinFit(const Matrix& fx, const Vector& y, Vector& c)
    2422{
    2523  int n = y.NElts();
    … …  
    2826  int nf = fx.NRows();
    2927
    30   OMatrix a(nf,nf);
     28  Matrix a(nf,nf);
    3129
    3230  for (int j=0; j<nf; j++)
    … …  
    3634  c = fx * y;
    3735
    38   if (OMatrix::GausPiv(a,c) == 0) THROW(singMatxErr);
     36  if (Matrix::GausPiv(a,c) == 0) THROW(singMatxErr);
    3937
    4038  double xi2 = 0;
    … …  
    5250
    5351
    54 double LinFit(const OVector& x, const OVector& y, const OVector& errY2, int nf,
    55            double (*f)(int, double), OVector& c, OVector& errC)
     52double LinFit(const Vector& x, const Vector& y, const Vector& errY2, int nf,
     53           double (*f)(int, double), Vector& c, Vector& errC)
    5654{
    5755  int n = x.NElts();
    5856  if (n != y.NElts()) THROW(sizeMismatchErr);
    5957 
    60   OMatrix fx(nf, n);
     58  Matrix fx(nf, n);
    6159  for (int i=0; i<nf; i++)
    6260    for (int j=0; j<n; j++)
    … …  
    6765
    6866
    69 double LinFit(const OMatrix& fx, const OVector& y, const OVector& errY2,
    70            OVector& c, OVector& errC)
     67double LinFit(const Matrix& fx, const Vector& y, const Vector& errY2,
     68           Vector& c, Vector& errC)
    7169{
    7270  int n = y.NElts();
    … …  
    7674  int nf = fx.NRows();
    7775
    78   OMatrix a(nf,nf);
     76  Matrix a(nf,nf);
    7977
    8078  c.Realloc(nf);
    … …  
    8987    }
    9088 
    91   OMatrix d(nf,nf+1);
     89  Matrix d(nf,nf+1);
    9290  for (int k=0; k<nf; k++) {
    9391    double x=0;
    … …  
    9997  }
    10098
    101   if (OMatrix::GausPiv(a,d) == 0) THROW(singMatxErr);
     99  if (Matrix::GausPiv(a,d) == 0) THROW(singMatxErr);
    102100
    103101  for (int l=0; l<nf; l++) {

  • trunk/SophyaLib/NTools/linfit.h

    r508 r514  
    11// This may look like C code, but it is really -*- C++ -*-
    22//
    3 // $Id: linfit.h,v 1.3 1999-10-25 10:36:08 ansari Exp $
     3// $Id: linfit.h,v 1.4 1999-10-25 16:40:00 ansari Exp $
    44//
    55
    … …  
    1010
    1111#include "machdefs.h"
    12 class OMatrix;
    13 class OVector;
     12#include "tvector.h"
    1413
    15 double LinFit(const OVector& x, const OVector& y, int nf,
    16            double (*f)(int, double), OVector& c);
     14namespace PlanckDPC {
     15
     16double LinFit(const Vector& x, const Vector& y, int nf,
     17           double (*f)(int, double), Vector& c);
    1718// fit lineaire des y en tant que somme de c(i)f(i,x), i=0..nf-1;
    1819
    1920
    20 double LinFit(const OMatrix& fx, const OVector& y, OVector& c);
     21double LinFit(const Matrix& fx, const Vector& y, Vector& c);
    2122// fit lineaire des y en tant que somme de c(i)f(i,x), i=0..nf-1,
    2223// la matrice fx contient les valeurs des f:
    … …  
    2425                     
    2526
    26 double LinFit(const OVector& x, const OVector& y, const OVector& errY2, int nf,
    27            double (*f)(int, double), OVector& c, OVector& errC);
     27double LinFit(const Vector& x, const Vector& y, const Vector& errY2, int nf,
     28           double (*f)(int, double), Vector& c, Vector& errC);
    2829// fit lineaire des y en tant que somme de c(i)f(i,x), i=0..nf-1,
    2930// errY2 contient les carres des erreurs sur les Y.
    … …  
    3132
    3233
    33 double LinFit(const OMatrix& fx, const OVector& y, const OVector& errY2,
    34            OVector& c, OVector& errC);
     34double LinFit(const Matrix& fx, const Vector& y, const Vector& errY2,
     35           Vector& c, Vector& errC);
    3536// fit lineaire des y en tant que somme de c(i)f(i,x), i=0..nf-1,
    3637// la matrice fx contient les valeurs des f:
    … …  
    3940// au retour, errC contient les erreurs sur les coefs.
    4041
     42} // Fin du namespace
     43
    4144#endif // LINFIT_SEEN

  • trunk/SophyaLib/NTools/outilsinit.cc

    r508 r514  
    3434  PPRegister(OMatrix);
    3535  PPRegister(OVector);
    36   PPRegister(Poly);
    37   PPRegister(Poly2);
     36  PPRegister(ObjFileIO<Poly>);
     37  PPRegister(ObjFileIO<Poly2>);
    3838  PPRegister(ObjFileIO<DVList>);
    3939  PPRegister(ObjFileIO<Histo>);

  • trunk/SophyaLib/NTools/poly.cc

    r508 r514  
    1313//++
    1414// Links        Parents
    15 // OVector
     15// Vector
    1616//--
    1717
    … …  
    2121
    2222
     23//////////////////////////////////////////////////////////////////////////
    2324//++
    2425// Poly::Poly(int degre=0)
    … …  
    2930
    3031Poly::Poly(int degre)
    31 : OVector(degre+1), dirty(0), deg(0)
     32: Vector(degre+1), dirty(0), deg(0)
    3233{
    3334  END_CONSTRUCTOR
    3435}
    3536
     37//++
     38Poly::Poly(Poly const& a)
     39//
     40//      Constructeur par copie.
     41//--
     42:Vector(a), dirty(a.dirty), deg(a.deg)
     43{
     44  END_CONSTRUCTOR
     45}
     46
    3647
    3748void Poly::UpdateDeg() const
    3849{
    39   int i = nalloc-1;
    40   while (data[i] == 0 && i>0) i--;
     50  int i = NElts()-1;
     51  while (Element(i) == 0 && i>0) i--;
    4152
    4253  (int&) deg = i;          // bientot mutable dans ANSI C++
    … …  
    6172{
    6273  UpdateDegIfDirty();
    63   double res = data[deg];
     74  double res = Element(deg);
    6475  for (int i=deg-1; i>=0; i--) {
    6576    res *= x;
    66     res += data[i];
     77    res += Element(i);
    6778  }
    6879  return res;
    … …  
    7687{
    7788  UpdateDegIfDirty();
    78   if (deg == 0) { data[0] = 0; return;}
     89  if (deg == 0) { Element(0) = 0; return;}
    7990  for (int i=1; i<=deg; i++)
    80     data[i-1] = data[i]*i;
    81   data[deg] = 0;
     91    Element(i-1) = Element(i)*i;
     92  Element(deg) = 0;
    8293  deg--;
    8394}
    … …  
    94105//  der.Zero();    // on sait que Realloc met a zero le reste.
    95106  for (int i=1; i<=deg; i++)
    96     der.data[i-1] = data[i]*i;
    97 }
    98 
    99 
    100 //++
    101 int Poly::Roots(OVector& roots) const
     107    der.Element(i-1) = Element(i)*i;
     108}
     109
     110
     111//++
     112int Poly::Roots(Vector& roots) const
    102113//
    103114//      Retourne dans roots les racines réelles, si on sait
    … …  
    133144  if (deg != 1) THROW(sizeMismatchErr);
    134145
    135   if (data[1] == 0) return 0;
    136   r = -data[0]/data[1];
     146  if (Element(1) == 0) return 0;
     147  r = -Element(0)/Element(1);
    137148  return 1;
    138149}
    … …  
    149160  if (deg != 2) THROW(sizeMismatchErr);
    150161
    151   double delta = data[1]*data[1] - 4*data[0]*data[2];
     162  double delta = Element(1)*Element(1) - 4*Element(0)*Element(2);
    152163  if (delta < 0) return 0;
    153164  if (delta == 0) {
    154     r1 = r2 = -data[1]/2/data[0];
     165    r1 = r2 = -Element(1)/2/Element(0);
    155166    return 1;
    156167  }
    157   r1 = (-data[1] - sqrt(delta))/2/data[0];
    158   r2 = (-data[1] + sqrt(delta))/2/data[0];
     168  r1 = (-Element(1) - sqrt(delta))/2/Element(0);
     169  r2 = (-Element(1) + sqrt(delta))/2/Element(0);
    159170  return 2;
    160171}
    … …  
    167178{
    168179  if (this == &a) return *this;
    169   OVector::operator=(a);
     180  Vector::operator=(a);
    170181
    171182  UpdateDeg();
    … …  
    185196  b.UpdateDegIfDirty();
    186197
    187   if (b.deg > deg)
    188     Realloc(b.deg);
     198  if (b.deg > deg) Realloc(b.deg);
    189199
    190200  int n = (deg > b.deg) ? deg : b.deg;
    191   for (int i=0; i<=n; i++)
    192     data[i] += b.data[i];
     201  for (int i=0; i<=n; i++) Element(i) += b.Element(i);
    193202
    194203  UpdateDeg();
    … …  
    204213  b.UpdateDegIfDirty();
    205214
    206   if (b.deg > deg)
    207     Realloc(b.deg);
     215  if (b.deg > deg) Realloc(b.deg);
    208216
    209217  int n = (deg > b.deg) ? deg : b.deg;
    210   for (int i=0; i<=n; i++)
    211     data[i] -= b.data[i];
     218  for (int i=0; i<=n; i++) Element(i) -= b.Element(i);
    212219
    213220  UpdateDeg();
    … …  
    215222}
    216223
    217 
    218 //++
    219 Poly operator+ (Poly const& a, Poly const& b)
    220 //
    221 //--
    222 #if HAS_NAMED_RETURN
    223      return c((a.deg > b.deg)?(a.deg+1):(b.deg+1))
    224 #endif
    225 {
    226 #if !HAS_NAMED_RETURN
    227 Poly c((a.deg > b.deg)?(a.deg+1):(b.deg+1));
    228 #endif
    229   c = a;
    230   c += b;
    231 #if !HAS_NAMED_RETURN
    232 return c;
    233 #endif
    234 }
    235 
    236 //++
    237 Poly operator- (Poly const& a, Poly const& b)
    238 //
    239 //--
    240 #if HAS_NAMED_RETURN
    241      return c((a.deg > b.deg)?(a.deg+1):(b.deg+1))
    242 #endif
    243 {
    244 #if !HAS_NAMED_RETURN
    245 Poly c((a.deg > b.deg)?(a.deg+1):(b.deg+1));
    246 #endif
    247   c = a;
    248   c -= b;
    249 #if !HAS_NAMED_RETURN
    250 return c;
    251 #endif
    252 }
    253 
    254 //++
    255 Poly operator* (Poly const& a, Poly const& b)
    256 //
    257 //--
    258 //#if HAS_NAMED_RETURN
    259 //     return c(a.deg + b.deg)
    260 //#endif
    261 {
    262 //#if !HAS_NAMED_RETURN
    263 Poly c(a.deg + b.deg);
    264 //#endif
    265   a.UpdateDegIfDirty();
     224//++
     225Poly& Poly::operator *= (double a)
     226//
     227//--
     228{
     229  UpdateDegIfDirty();
     230  for (int i=0; i<=deg; i++) Element(i) *= a;
     231  return *this;
     232}
     233
     234//++
     235Poly Poly::Mult(Poly const& b) const
     236//
     237//--
     238{
     239Poly c(deg + b.deg);
     240  UpdateDegIfDirty();
    266241  b.UpdateDegIfDirty();
    267242
    268   c.deg = a.deg + b.deg;
     243  c.deg = deg + b.deg;
    269244
    270245  for (int i=0; i<=c.deg; i++) {
    271246    c[i] = 0;
    272     int kmin = (i <= a.deg) ? 0 : i - a.deg;
    273     int kmax = (i <= a.deg) ? i : a.deg;
     247    int kmin = (i <= deg) ? 0 : i - deg;
     248    int kmax = (i <= deg) ? i : deg;
    274249    for (int k=kmin; k<=kmax; k++)
    275       c[i] += a[k] * b[i-k];
    276   }
    277 //#if !HAS_NAMED_RETURN
     250      c[i] += (*this)[k] * b[i-k];
     251  }
    278252return c;
    279 //#endif
    280 }
    281 
    282 //++
    283 Poly& Poly::operator *= (double a)
    284 //
    285 //--
    286 {
    287   UpdateDegIfDirty();
    288 
    289   for (int i=0; i<=deg; i++)
    290     data[i] *= a;
    291 
    292   return *this;
    293 }
    294 
    295 
    296 //++
    297 Poly operator* (double a, Poly const& b)
    298 //
    299 //--
    300 #if HAS_NAMED_RETURN
    301      return c(b)
    302 #endif
    303 {
    304 #if !HAS_NAMED_RETURN
    305   Poly c(b);
    306 #endif
    307 
    308    c *= a;
    309 
    310 #if !HAS_NAMED_RETURN
    311   return c;
    312 #endif
    313 }
    314 
    315 
    316 void Poly::ReadSelf(PInPersist& s)
    317 {
    318   int_4 dg;
    319   s >> dg;
    320   Realloc(dg, true);
    321 
    322 #ifdef DEBUG
    323   int r = nr; int c = nc;
    324 #endif
    325 
    326   OVector::ReadSelf(s);
    327 
    328 #ifdef DEBUG
    329   DBASSERT(r == nr && c == nc);
    330 #endif
    331   UpdateDeg();
    332 
    333   #ifdef DEBUG
    334   DBASSERT(dg == deg);
    335   #endif
    336 }
    337 
    338 void Poly::WriteSelf(POutPersist& s) const
    339 {
    340   UpdateDegIfDirty();
    341   ((Poly*)(this))->Realloc(deg, true);
    342   s << deg;
    343   OVector::WriteSelf(s);
    344253}
    345254
    … …  
    365274
    366275//++
    367 ostream& operator << (ostream& s, const Poly& a)
    368 //
    369 //      Impressions sur le stream "s".
    370 //--
    371 {
    372   a.Print(s);
    373   return s;
    374 }
    375 
    376 //++
    377 double Poly::Fit(OVector const& x, OVector const& y, int degre)
     276double Poly::Fit(Vector const& x, Vector const& y, int degre)
    378277//
    379278//      Ajustement polynomial par moindre carrés. Un polynôme de
    … …  
    387286  Realloc(degre);
    388287
    389   OMatrix a(degre+1, n);
     288  Matrix a(degre+1, n);
    390289
    391290  for (int c=0; c<n; c++) {
    … …  
    397296  }
    398297
    399   double rc = LinFit(a,y,(OVector&)*this);
     298  double rc = LinFit(a,y,(Vector&)*this);
    400299  UpdateDeg();
    401300  return rc;
    … …  
    403302
    404303//++
    405 double Poly::Fit(OVector const& x, OVector const& y, OVector const& erry2, int degre,
    406                  OVector& errCoef)
     304double Poly::Fit(Vector const& x, Vector const& y, Vector const& erry2, int degre,
     305                 Vector& errCoef)
    407306//
    408307//      Ajustement polynomial par moindre carrés. Un polynôme de
    … …  
    419318  errCoef.Realloc(degre+1);
    420319
    421   OMatrix a(degre+1, n);
     320  Matrix a(degre+1, n);
    422321
    423322  for (int c=0; c<n; c++) {
    … …  
    429328  }
    430329
    431   double rc = LinFit(a,y,erry2,(OVector&)*this,errCoef);
     330  double rc = LinFit(a,y,erry2,(Vector&)*this,errCoef);
    432331  UpdateDeg();
    433332  return rc;
    434333}
    435334
     335
     336//++
     337// Poly Poly::power(int n) const
     338//
     339//      Retourne le polynôme à la puissance n
     340//--
     341
     342Poly Poly::power(int n) const // a accelerer !!!
     343{
     344  if (n < 0) THROW(rangeCheckErr);
     345  if (n == 0) { Poly r(0); r[0] = 1; return r;}
     346  if (n == 1) { return *this; }
     347  return *this * power(n-1);
     348}
     349
     350//++
     351Poly Poly::operator() (Poly const& b) const
     352//
     353//      Substitution d'un polynôme dans un autre.
     354//--
     355{
     356  Poly c(b.Degre()*Degre());
     357  for (int i=0; i<= Degre(); i++)
     358    c += (*this)[i] * b.power(i);
     359  return c;
     360}
     361
     362
     363//////////////////////////////////////////////////////////////////////////
     364void ObjFileIO<Poly>::ReadSelf(PInPersist& is)
     365{
     366if(dobj==NULL) dobj=new Poly;
     367int_4 dg;
     368is >> dg;
     369dobj->Realloc(dg,true);
     370is >> *((Vector *) dobj);
     371dobj->UpdateDeg();
     372}
     373
     374void ObjFileIO<Poly>::WriteSelf(POutPersist& os) const
     375{
     376if(dobj == NULL) return;
     377dobj->UpdateDegIfDirty();
     378dobj->Realloc(dobj->deg,true);
     379os << dobj->deg;
     380os << *((Vector *) dobj);
     381}
     382
     383//////////////////////////////////////////////////////////////////////////
    436384int binomial(int n, int p)
    437385{
    … …  
    442390}
    443391
    444 
    445 //++
    446 // Poly Poly::power(int n) const
    447 //
    448 //      Retourne le polynôme à la puissance n
    449 //--
    450 
    451 Poly Poly::power(int n) const // a accelerer !!!
    452 {
    453   if (n < 0) THROW(rangeCheckErr);
    454   if (n == 0) { Poly r(0); r[0] = 1; return r;}
    455   if (n == 1) { return *this; }
    456   return *this * power(n-1);
    457 }
    458 
    459 //++
    460 Poly Poly::operator() (Poly const& b) const
    461 //
    462 //      Substitution d'un polynôme dans un autre.
    463 //--
    464 {
    465   Poly c(b.Degre()*Degre());
    466   for (int i=0; i<= Degre(); i++)
    467     c += (*this)[i] * b.power(i);
    468   return c;
    469 }
    470 
    471 
     392//////////////////////////////////////////////////////////////////////////
    472393// ******************* POLY 2 VARIABLES ******************
    473 
    474394//++
    475395// Class        Poly2
    … …  
    482402//++
    483403// Links        Parents
    484 // OVector
     404// Vector
    485405//--
    486406
    … …  
    494414//      Crée un polynôme de degrés partiels degreX et degreY.
    495415//--
    496 :OVector((degreX+1)*(degreY+1)), dirty(0),
     416:Vector((degreX+1)*(degreY+1)), dirty(0),
    497417 maxDegX(degreX), maxDegY(degreY), degX(0), degY(0), deg(0)
    498418{
    … …  
    506426//      de deux polynômes à une variable, p2(x,y)=px(x)py(y)
    507427//--
    508 :OVector((polX.Degre()+1)*(polY.Degre()+1)), dirty(0),
     428:Vector((polX.Degre()+1)*(polY.Degre()+1)), dirty(0),
    509429 maxDegX(polX.Degre()), maxDegY(polY.Degre()),
    510430 degX(polX.Degre()), degY(polY.Degre()), deg(degX+degY)
    … …  
    521441//      Constructeur par copie.
    522442//--
    523 :OVector(a), dirty(a.dirty),
     443:Vector(a), dirty(a.dirty),
    524444 maxDegX(a.maxDegX), maxDegY(a.maxDegY),
    525445 degX(a.degX), degY(a.degY), deg(a.deg)
    … …  
    561481  UpdateDegIfDirty();
    562482  Poly2 tmp(*this);
    563   OVector::Realloc((degreX+1)*(degreY+1));
    564   Zero();
     483  Vector::Realloc((degreX+1)*(degreY+1));
     484  Reset();
    565485  maxDegX = degreX;
    566486  maxDegY = degreY;
    … …  
    631551
    632552//++
    633 double Poly2::Fit(OVector const& x, OVector const& y, OVector const& z,
     553double Poly2::Fit(Vector const& x, Vector const& y, Vector const& z,
    634554                  int degreX, int degreY)
    635555//
    … …  
    643563  Realloc(degreX, degreY);
    644564
    645   OMatrix a((degreX+1)*(degreY+1), n);
     565  Matrix a((degreX+1)*(degreY+1), n);
    646566
    647567  for (int c=0; c<n; c++) {
    … …  
    657577  }
    658578
    659   double rc = LinFit(a,z,(OVector&)*this);
     579  double rc = LinFit(a,z,(Vector&)*this);
    660580  UpdateDeg();
    661581  return rc;
    … …  
    664584
    665585//++
    666 double Poly2::Fit(OVector const& x, OVector const& y, OVector const& z,
    667                   OVector const& errz2, int degreX, int degreY,
    668                   OVector& errCoef)
     586double Poly2::Fit(Vector const& x, Vector const& y, Vector const& z,
     587                  Vector const& errz2, int degreX, int degreY,
     588                  Vector& errCoef)
    669589//
    670590//      Ajustement par moindre carrés z = P(x,y), degrés partiels imposés,
    … …  
    680600  errCoef.Realloc((degreX+1)*(degreY+1));
    681601
    682   OMatrix a((degreX+1)*(degreY+1), n);
     602  Matrix a((degreX+1)*(degreY+1), n);
    683603
    684604  for (int c=0; c<n; c++) {
    … …  
    694614  }
    695615
    696   double rc = LinFit(a,z,errz2,(OVector&)*this,errCoef);
     616  double rc = LinFit(a,z,errz2,(Vector&)*this,errCoef);
    697617  UpdateDeg();
    698618  return rc;
    … …  
    700620
    701621//++
    702 double Poly2::Fit(OVector const& x, OVector const& y, OVector const& z,
     622double Poly2::Fit(Vector const& x, Vector const& y, Vector const& z,
    703623                  int degre)
    704624//
    … …  
    712632  Realloc(degre, degre);   // certains vaudront 0, impose.
    713633
    714   OMatrix a((degre+1)*(degre+2)/2, n);
    715   OVector cf((degre+1)*(degre+2)/2);
     634  Matrix a((degre+1)*(degre+2)/2, n);
     635  Vector cf((degre+1)*(degre+2)/2);
    716636#define C_INDEX(i,j) ((i) + (j)*(2*degre+3-(j))/2)
    717637
    … …  
    744664
    745665//++
    746 double Poly2::Fit(OVector const& x, OVector const& y, OVector const& z,
    747                   OVector const& errz2, int degre,
    748                   OVector& errCoef)
     666double Poly2::Fit(Vector const& x, Vector const& y, Vector const& z,
     667                  Vector const& errz2, int degre,
     668                  Vector& errCoef)
    749669//
    750670//      Ajustement par moindre carrés z = P(x,y), degré total imposé,
    … …  
    761681#define C_INDEX(i,j) ((i) + (j)*(2*degre+3-(j))/2)
    762682
    763   OMatrix a((degre+1)*(degre+2)/2, n);
    764   OVector cf((degre+1)*(degre+2)/2);
    765   OVector ecf((degre+1)*(degre+2)/2);
     683  Matrix a((degre+1)*(degre+2)/2, n);
     684  Vector cf((degre+1)*(degre+2)/2);
     685  Vector ecf((degre+1)*(degre+2)/2);
    766686
    767687  for (int c=0; c<n; c++) {
    … …  
    792712  UpdateDeg();
    793713  return rc;
    794 }
    795 
    796 
    797 void Poly2::ReadSelf(PInPersist& s)
    798 {
    799   int_4 dgx, dgy;
    800   s >> dgx >> dgy;
    801 
    802   Realloc(dgx, dgy);
    803 
    804 #ifdef DEBUG
    805   int r = nr; int c = nc;
    806 #endif
    807 
    808   OVector::ReadSelf(s);
    809 
    810 #ifdef DEBUG
    811   DBASSERT(r == nr && c == nc);
    812 #endif
    813   UpdateDeg();
    814 }
    815 
    816 void Poly2::WriteSelf(POutPersist& s) const
    817 {
    818   s << maxDegX << maxDegY;
    819   OVector::WriteSelf(s);
    820714}
    821715
    … …  
    850744
    851745//++
    852 ostream& operator << (ostream& s, const Poly2& a)
    853 //
    854 //      Impression sur stream s.
    855 //--
    856 {
    857   a.Print(s);
    858   return s;
    859 }
    860 
    861 
    862 //++
    863746// Titre        Opérations
    864747//--
    … …  
    905788
    906789//++
    907 Poly2 operator+ (Poly2 const& a, Poly2 const& b)
    908 //
    909 //--
    910 #if HAS_NAMED_RETURN
    911      return c(a)
    912 #endif
    913 {
    914 #if !HAS_NAMED_RETURN
    915 Poly2 c(a);
    916 #endif
    917   c += b;
    918 #if !HAS_NAMED_RETURN
    919 return c;
    920 #endif
    921 }
    922 
    923 //++
    924 Poly2 operator- (Poly2 const& a, Poly2 const& b)
    925 //
    926 //--
    927 #if HAS_NAMED_RETURN
    928      return c(a)
    929 #endif
    930 {
    931 #if !HAS_NAMED_RETURN
    932 Poly2 c(a);
    933 #endif
    934   c -= b;
    935 #if !HAS_NAMED_RETURN
    936 return c;
    937 #endif
    938 }
    939 
    940 //++
    941790Poly2& Poly2::operator *= (double a)
    942791//
    … …  
    944793{
    945794  for (int i=0; i<NElts(); i++)
    946     data[i] *= a;
     795    Element(i) *= a;
    947796
    948797  return *this;
    … …  
    950799
    951800//++
    952 Poly2 operator * (double a, Poly2 const& b)
    953 //
    954 //--
    955 #if HAS_NAMED_RETURN
    956      return c(b)
    957 #endif
    958 {
    959 #if !HAS_NAMED_RETURN
    960 Poly2 c(b);
    961 #endif
    962   c *= a;
    963 #if !HAS_NAMED_RETURN
    964 return c;
    965 #endif
    966 }
    967 
    968 //++
    969 Poly2 operator* (Poly2 const& a, Poly2 const& b)
    970 //
    971 //--
    972 {
    973   Poly2 c(a.DegX() + b.DegX(), a.DegY() + b.DegY());
    974   a.UpdateDegIfDirty();
     801Poly2 Poly2::Mult(Poly2 const& b) const
     802//
     803//--
     804{
     805  Poly2 c(DegX() + b.DegX(), DegY() + b.DegY());
     806  UpdateDegIfDirty();
    975807  b.UpdateDegIfDirty();
    976808
    977   for (int i=0; i<=a.DegX(); i++)
    978     for (int j=0; j<=a.DegY(); j++)
     809  for (int i=0; i<=DegX(); i++)
     810    for (int j=0; j<=DegY(); j++)
    979811      for (int k=0; k<=b.DegX(); k++)
    980812        for (int l=0; l<=b.DegY(); l++)
    981           c.Coef(i+k,j+l) += a.Coef(i,j)*b.Coef(k,l);
     813          c.Coef(i+k,j+l) += Coef(i,j)*b.Coef(k,l);
    982814  return c;
    983815}
    … …  
    1029861}
    1030862
    1031 
     863//////////////////////////////////////////////////////////////////////////
     864void ObjFileIO<Poly2>::ReadSelf(PInPersist& is)
     865{
     866if(dobj==NULL) dobj=new Poly2;
     867int_4 dgx, dgy;
     868is >> dgx >> dgy;
     869dobj->Realloc(dgx,dgy);
     870is >> *((Vector *) dobj);
     871dobj->UpdateDeg();
     872}
     873
     874void ObjFileIO<Poly2>::WriteSelf(POutPersist& os) const
     875{
     876if(dobj == NULL) return;
     877os << dobj->maxDegX << dobj->maxDegY;
     878os << *((Vector *) dobj);
     879}
     880
     881
     882//////////////////////////////////////////////////////////////////////////
     883#ifdef __CXX_PRAGMA_TEMPLATES__
     884#pragma define_template ObjFileIO<Poly>
     885#pragma define_template ObjFileIO<Poly2>
     886#endif
     887
     888#if defined(ANSI_TEMPLATES) || defined(GNU_TEMPLATES)
     889template class ObjFileIO<Poly>;
     890template class ObjFileIO<Poly2>;
     891#endif

  • trunk/SophyaLib/NTools/poly.h

    r508 r514  
    11// This may look like C code, but it is really -*- C++ -*-
    22//
    3 // $Id: poly.h,v 1.2 1999-10-25 10:36:11 ansari Exp $
     3// $Id: poly.h,v 1.3 1999-10-25 16:40:01 ansari Exp $
    44//
    5 
    65
    76// Des polynomes avec des operations de bases et surtout des fits.
    … …  
    109#define POLY_SEEN
    1110
     11#include "objfio.h"
    1212#include <stdio.h>
    1313#include "peida.h"
    14 #include "cvector.h"
     14#include "tvector.h"
     15#include "ppersist.h"
     16#include "anydataobj.h"
     17
     18namespace PlanckDPC {
    1519
    1620class Poly2;
    1721
    18 class Poly : public OVector {
     22//////////////////////////////////////////////////////////////////////////
     23class Poly : public Vector {
     24  friend class ObjFileIO<Poly>;
    1925public:
    2026  Poly(int degre = 0);
     27  Poly(Poly const& a);
    2128
    2229  enum {classId = ClassId_Poly1};
    23   int_4                 ClassId() const         { return classId; }
    24   static PPersist*      Create() {return new Poly;}
     30  int_4  ClassId() const         { return classId; }
    2531
    2632  inline int Degre() const {UpdateDegIfDirty(); return deg;}
    2733
    28   inline void Realloc(int n, bool force=false) {OVector::Realloc(n+1,force);}
    29 
    30   inline double operator[](int i) const {return data[i];}
    31   inline double& operator[](int i) {dirty = 1; return data[i];}
     34  inline void Realloc(int n, bool force=false) {Vector::Realloc(n+1,force);}
     35
     36  inline double operator[](int i) const {return Element(i);}
     37  inline double& operator[](int i) {dirty = 1; return Element(i);}
    3238  // Retourne le coefficient de degre i
    3339
    … …  
    4147  // Derive le polynome dans un autre
    4248
    43   int    Roots(OVector& roots) const;
     49  int    Roots(Vector& roots) const;
    4450  // retourne les racines si on peut les calculer...
    4551
    … …  
    4955  int    Root2(double& r1, double& r2) const;
    5056  // special degre 2
    51 
    52   friend Poly operator* (Poly const& a, Poly const& b);
    53   friend Poly operator+ (Poly const& a, Poly const& b);
    54   friend Poly operator- (Poly const& a, Poly const& b);
    5557
    5658  Poly& operator = (Poly const& b);
    … …  
    5860  Poly& operator -= (Poly const& b);
    5961
    60   friend Poly operator* (double a, Poly const& b);
    61 
    6262  Poly& operator *= (double a);
     63
     64  Poly Mult(Poly const& b) const;
    6365
    6466  Poly power(int n) const;
    … …  
    6668  Poly2 operator() (Poly2 const& b) const;
    6769
    68   void ReadSelf(PInPersist&);
    69   void WriteSelf(POutPersist&) const;
    70 
    7170  void Print(ostream& s) const;
    72   friend ostream& operator << (ostream& s, const Poly& a);
    73 
    74   double Fit(OVector const& x, OVector const& y, int degre);
     71
     72  double Fit(Vector const& x, Vector const& y, int degre);
    7573  // Fit d'un polynome de degre donne sur les x et y.
    7674
    77   double Fit(OVector const& x, OVector const& y, OVector const& erry2, int degre,
    78             OVector& errCoef);
     75  double Fit(Vector const& x, Vector const& y, Vector const& erry2, int degre,
     76            Vector& errCoef);
    7977  // En plus, on fournit les carres des erreurs sur y et on a les erreurs
    8078  // sur les coefficients dans un vecteur.
    … …  
    8785};
    8886
     87inline Poly operator* (Poly const& a, Poly const& b)
     88  { return a.Mult(b); }
     89
     90inline Poly operator+ (Poly const& a, Poly const& b)
     91 { Poly c((a.Degre() > b.Degre())?(a.Degre()+1):(b.Degre()+1));
     92   c = a; c += b; return c; }
     93
     94inline Poly operator- (Poly const& a, Poly const& b)
     95  { Poly c((a.Degre() > b.Degre())?(a.Degre()+1):(b.Degre()+1));
     96    c = a; c -= b; return c; }
     97
     98inline Poly operator* (double a, Poly const& b)
     99  { Poly c(b); c *= a; return c; }
     100
     101inline ostream& operator << (ostream& s, const Poly& a)
     102  { a.Print(s); return s; }
     103
     104//////////////////////////////////////////////////////////////////////////
     105inline POutPersist& operator << (POutPersist& os, Poly & obj)
     106  { ObjFileIO<Poly> fio(&obj);  fio.Write(os);  return(os); }
     107inline PInPersist& operator >> (PInPersist& is, Poly & obj)
     108  { ObjFileIO<Poly> fio(&obj);  fio.Read(is);  return(is); }
     109// Classe pour la gestion de persistance
     110// ObjFileIO<Poly>
     111
     112//////////////////////////////////////////////////////////////////////////
    89113int binomial(int n, int p);
    90114
    91 
    92 class Poly2 : public OVector {  // Ca pourrait etre une matrice mais
    93                                   // la encore, un vecteur est utile pour les
    94                                   // fits.
    95 
     115//////////////////////////////////////////////////////////////////////////
     116class Poly2 : public Vector {
     117  friend class ObjFileIO<Poly2>;
    96118public:
    97119  Poly2(int degreX=0, int degreY=0);
    … …  
    101123
    102124  enum {classId = ClassId_Poly2};
    103   int_4                 ClassId() const         { return classId; }
    104   static PPersist*      Create() {return new Poly2;}
     125  int_4 ClassId() const         { return classId; }
    105126
    106127  inline int DegX() const {UpdateDegIfDirty(); return degX;}
    … …  
    123144
    124145  inline double Coef(int dx, int dy) const {
    125     return (dx>maxDegX || dy>maxDegY) ? 0 : data[IndCoef(dx,dy)];
     146    return (dx>maxDegX || dy>maxDegY) ? 0 : Element(IndCoef(dx,dy));
    126147  }
    127148  inline double& Coef(int dx, int dy) {
    128149    if (dx>maxDegX || dy>maxDegY) THROW(rangeCheckErr);
    129     dirty = 1; return data[IndCoef(dx,dy)];
     150    dirty = 1; return Element(IndCoef(dx,dy));
    130151  }
    131152  // retourne le coefficient de degre (dx,dy)
    132153
    133   double Fit(OVector const& x, OVector const& y, OVector const& z,
     154  double Fit(Vector const& x, Vector const& y, Vector const& z,
    134155             int degreX, int degreY);
    135   double Fit(OVector const& x, OVector const& y, OVector const& z,
    136              OVector const& errz2, int degreX, int degreY,
    137              OVector& errCoef);
     156  double Fit(Vector const& x, Vector const& y, Vector const& z,
     157             Vector const& errz2, int degreX, int degreY,
     158             Vector& errCoef);
    138159  // degres partiels imposes. cf Poly::Fit sinon
    139160
    140161
    141   double Fit(OVector const& x, OVector const& y, OVector const& z,
     162  double Fit(Vector const& x, Vector const& y, Vector const& z,
    142163             int degre);
    143   double Fit(OVector const& x, OVector const& y, OVector const& z,
    144              OVector const& errz2, int degre,
    145              OVector& errCoef);
     164  double Fit(Vector const& x, Vector const& y, Vector const& z,
     165             Vector const& errz2, int degre,
     166             Vector& errCoef);
    146167  // degre total impose. cf Poly::Fit sinon
    147168
    … …  
    151172  Poly2& operator -= (Poly2 const& b);
    152173
    153   friend Poly2 operator + (Poly2 const& a, Poly2 const& b);
    154   friend Poly2 operator - (Poly2 const& a, Poly2 const& b);
    155   friend Poly2 operator * (Poly2 const& a, Poly2 const& b);
     174  Poly2& operator *= (double a);
     175
     176  Poly2 Mult(Poly2 const& b) const;
    156177
    157178  Poly2 power(int n) const;
    158 
    159   Poly2& operator *= (double a);
    160   friend Poly2 operator * (double a, Poly2 const& b);
    161179
    162180  Poly2 operator() (Poly const& px, Poly const& py) const;
    … …  
    165183  void Realloc(int degreX, int degreY);
    166184
    167   void ReadSelf(PInPersist&);
    168   void WriteSelf(POutPersist&) const;
    169 
    170185  void Print(ostream& s) const;
    171   friend ostream& operator << (ostream& s, const Poly2& a);
    172186
    173187private:
    … …  
    182196};
    183197
     198inline Poly2 operator* (Poly2 const& a, Poly2 const& b)
     199  { return a.Mult(b); }
     200
     201inline Poly2 operator+ (Poly2 const& a, Poly2 const& b)
     202  { Poly2 c(a); c += b; return c; }
     203
     204inline Poly2 operator- (Poly2 const& a, Poly2 const& b)
     205  { Poly2 c(a); c -= b; return c; }
     206
     207inline Poly2 operator * (double a, Poly2 const& b)
     208  { Poly2 c(b); c *= a; return c; }
     209
     210inline ostream& operator << (ostream& s, const Poly2& a)
     211  { a.Print(s); return s; }
     212
     213//////////////////////////////////////////////////////////////////////////
     214inline POutPersist& operator << (POutPersist& os, Poly2 & obj)
     215  { ObjFileIO<Poly2> fio(&obj);  fio.Write(os);  return(os); }
     216inline PInPersist& operator >> (PInPersist& is, Poly2 & obj)
     217  { ObjFileIO<Poly2> fio(&obj);  fio.Read(is);  return(is); }
     218// Classe pour la gestion de persistance
     219// ObjFileIO<Poly2>
     220
     221
     222} // Fin du namespace
     223
    184224#endif // POLY_SEEN

  • trunk/SophyaLib/NTools/rk4cdifeq.cc

    r508 r514  
    7474//++
    7575RK4CDiffEq& 
    76 RK4CDiffEq::AbsAccuracy(OVector const& vScal)
     76RK4CDiffEq::AbsAccuracy(Vector const& vScal)
    7777//
    7878//      On souhaite une précision absolue, et le vecteur vScal contient
    … …  
    118118
    119119void
    120 RK4CDiffEq::RKCStep(OVector& newY, OVector const& y0, OVector const& yScale,
     120RK4CDiffEq::RKCStep(Vector& newY, Vector const& y0, Vector const& yScale,
    121121               double dttry, double& dtdone, double& dtnext)
    122122{
    … …  
    151151  // Et on corrige a l'ordre 5
    152152
    153   newY += yTemp/15;
     153  newY += yTemp/15.;
    154154}
    155155
    156156void
    157 RK4CDiffEq::SolveArr(OMatrix& y, double* t, double tf, int n)
     157RK4CDiffEq::SolveArr(Matrix& y, double* t, double tf, int n)
    158158{
    159159  //TIMEF;
    … …  
    162162  double dxres = (tf - mXStart)/n;
    163163
    164   OVector yt = mYStart;
     164  Vector yt = mYStart;
    165165 
    166166  k1.Realloc(mFunc->NFuncReal());

  • trunk/SophyaLib/NTools/rk4cdifeq.h

    r508 r514  
    44
    55#include "difeq.h"
     6
     7namespace PlanckDPC {
    68
    79// <summary> Runge-Kutta ordre 4 adaptatif </summary>
    … …  
    2022  // <group>
    2123  RK4CDiffEq& Accuracy(double);
    22   RK4CDiffEq& AbsAccuracy(OVector const& vScal);
     24  RK4CDiffEq& AbsAccuracy(Vector const& vScal);
    2325  RK4CDiffEq& AbsAccuracy(double scal);
    2426  RK4CDiffEq& RelAccuracy();
    … …  
    2628 
    2729  // Implementation RK4 adaptatif
    28   virtual void SolveArr(OMatrix& y, double* t, double tf, int n);
     30  virtual void SolveArr(Matrix& y, double* t, double tf, int n);
    2931
    3032protected:
    3133  // Un pas adaptatif
    32   void RKCStep(OVector& newY, OVector const& y0, OVector const& yScale,
     34  void RKCStep(Vector& newY, Vector const& y0, Vector const& yScale,
    3335               double dttry, double& dtdone, double& dtnext);
    3436
    3537  double eps;
    3638  bool   relAccuracy;
    37   OVector accScale;
    38   OVector yTemp; // Pour ne pas reallouer
    39   OVector ySave;
     39  Vector accScale;
     40  Vector yTemp; // Pour ne pas reallouer
     41  Vector ySave;
    4042};
    4143
     44} // Fin du namespace
     45
    4246#endif

  • trunk/SophyaLib/NTools/rzimage.cc

    r508 r514  
    55// LAL (Orsay) / IN2P3-CNRS  DAPNIA/SPP (Saclay) / CEA
    66
    7 //  $Id: rzimage.cc,v 1.5 1999-10-25 10:36:12 ansari Exp $
     7//  $Id: rzimage.cc,v 1.6 1999-10-25 16:40:03 ansari Exp $
    88
    99#include "machdefs.h"
    … …  
    765765GeneralFunction* f = gfit.GetFunction();
    766766if(f==NULL) return NULL;
    767 OVector par = gfit.GetParm();
     767Vector par = gfit.GetParm();
    768768RzImage* img = new RzImage(DataType((r_4)0.),XSize(),YSize());
    769769img->SetOrg(XOrg(),YOrg());
    … …  
    790790GeneralFunction* f = gfit.GetFunction();
    791791if(f==NULL) return NULL;
    792 OVector par = gfit.GetParm();
     792Vector par = gfit.GetParm();
    793793RzImage* img = new RzImage(DataType((r_4)0.),XSize(),YSize());
    794794img->SetOrg(XOrg(),YOrg());

  • trunk/SophyaLib/NTools/tmatrix.cc

    r508 r514  
    1 // $Id: tmatrix.cc,v 1.11 1999-10-25 10:36:12 ansari Exp $
     1// $Id: tmatrix.cc,v 1.12 1999-10-25 16:40:03 ansari Exp $
    22//                         C.Magneville          04/99
    33#include "machdefs.h"
    … …  
    1010#include "tmatrix.h"
    1111#include "objfio.h"
    12 
    13 using namespace PlanckDPC;
    1412
    1513////////////////////////////////////////////////////////////////
    … …  
    184182}
    185183
    186 #include "matrix.h"
    187184////////////////////////////////////////////////////////////////
    188185//**** Pour inversion

  • trunk/SophyaLib/NTools/tmatrix.h

    r508 r514  
    5050  inline uint_4 NRows() const {return mNr;}
    5151  inline uint_4 NCols() const {return mNc;}
     52  inline uint_4 NCol() const {return mNc;} // back-compat Peida
    5253  inline T const& operator()(uint_4 r,uint_4 c) const
    5354                            {return *(mNDBlock.Begin()+r*mNc+c);}
    … …  
    139140
    140141////////////////////////////////////////////////////////////////
    141 // Surcharge d'operateurs A (+=,-=,*=,/=) (T) x
     142// Surcharge d'operateurs A (+,-,*,/) (T) x
    142143
    143144template <class T> inline TMatrix<T> operator + (const TMatrix<T>& a, T b)

  • trunk/SophyaLib/NTools/tvector.cc

    r306 r514  
    1 // $Id: tvector.cc,v 1.4 1999-05-19 10:17:42 ansari Exp $
     1// $Id: tvector.cc,v 1.5 1999-10-25 16:40:04 ansari Exp $
    22//                         C.Magneville          04/99
    33#include "machdefs.h"
    … …  
    99#include "tvector.h"
    1010#include "objfio.h"
    11 
    12 using namespace PlanckDPC;
    1311
    1412template <class T>

  • trunk/SophyaLib/NTools/tvector.h

    r508 r514  
    4747
    4848// produit scalaire, matrice*vecteur
    49 template <class T> inline T operator* (const TVector<T>& v1, const TVector<T>& v2)
     49template <class T>
     50inline T operator* (const TVector<T>& v1, const TVector<T>& v2)
    5051  {if(v1.NRows() != v2.NRows())
    5152     throw(SzMismatchError("TVector::operator*(TVector& v1,TVector v2) size mismatch"));
    … …  
    5556   return r;}
    5657
    57 template <class T> inline TVector<T> operator* (const TMatrix<T>& a, const TVector<T>& b)
    58                    {return TVector<T>(a * ((TMatrix<T> const&)(b)));}
     58template <class T>
     59inline TVector<T> operator* (const TMatrix<T>& a, const TVector<T>& b)
     60   {return TVector<T>(a * ((TMatrix<T> const&)(b)));}
    5961
    6062// Resolution du systeme A*C = B

  • trunk/SophyaLib/Samba/anagen.cc

    r508 r514  
    3737    double phi0;
    3838    float theta;
    39     OVector phin, datan, phis, datas;
     39    Vector phin, datan, phis, datas;
    4040    map.GetThetaSlice(map.NbThetaSlices()-ith-1,theta,phis,datas);
    4141    map.GetThetaSlice(ith,theta,phin,datan);
    … …  
    7777
    7878
    79 void comp_phas2gen(int nsmax,int nlmax,int nmmax, OVector datain,
     79void comp_phas2gen(int nsmax,int nlmax,int nmmax, Vector datain,
    8080                int nph,vector< complex<double> >& dataout,double phi0){
    8181  /*=======================================================================
    … …  
    9292  //FFTVector outf=FFTServer::solve(inf,-1);
    9393  FFTServer fft;
    94   OVector outf;
     94  Vector outf;
    9595  //cout<<"in :"<<datain<<endl;
    9696  fft.fftb(datain,outf);

  • trunk/SophyaLib/Samba/anagen.h

    r508 r514  
    44#include <complex>
    55#include <vector>
    6 #include "cvector.h"
     6#include "tvector.h"
    77#include "sphericalmap.h"
    88
    … …  
    1010             vector< vector< complex<double> > >& alm, double cos_theta_cut);
    1111
    12 void comp_phas2gen(int nsmax,int nlmax,int nmmax, OVector datain,
     12void comp_phas2gen(int nsmax,int nlmax,int nmmax, Vector datain,
    1313                int nph,vector< complex<double> >& dataout,double phi0);
    1414

  • trunk/SophyaLib/Samba/syngen.cc

    r508 r514  
    4343    double phi0;
    4444    float theta,thetas;
    45     OVector phin, datan,phis,datas;
     45    Vector phin, datan,phis,datas;
    4646    map.GetThetaSlice(map.NbThetaSlices()-ith-1,thetas,phis,datas);
    4747    map.GetThetaSlice(ith,theta,phin,datan);
    … …  
    104104void syn_phasgen(const int nsmax,const int nlmax,const int nmmax,
    105105              const vector< complex<long double> >& datain,
    106               int nph, OVector& dataout, double phi0,
     106              int nph, Vector& dataout, double phi0,
    107107              vector< complex<long double> >& bw){
    108108

  • trunk/SophyaLib/Samba/syngen.h

    r508 r514  
    44#include <complex>
    55#include <vector>
    6 #include "cvector.h"
     6#include "tvector.h"
    77#include "sphericalmap.h"
    88
    … …  
    1818void syn_phasgen(const int nsmax,const int nlmax,const int nmmax,
    1919              const vector< complex<long double> >& datain,
    20               int nph, OVector& dataout, double phi0,
     20              int nph, Vector& dataout, double phi0,
    2121              vector< complex<long double> >& bw);
    2222
Note: See TracChangeset for help on using the changeset viewer.