Changeset 3572 in Sophya for trunk/SophyaLib/NTools


Ignore:
Timestamp:
Feb 7, 2009, 10:50:34 PM (17 years ago)
Author:
cmv
Message:

char* -> const char* pour regler les problemes de deprecated string const... + comparaison unsigned signed + suppression EVOL_PLANCK rz+cmv 07/02/2009

Location:
trunk/SophyaLib/NTools
Files:
14 edited

Legend:

Unmodified
Added
Removed
  • trunk/SophyaLib/NTools/Makefile

    r3406 r3572  
    44 
    55clean:
    6          rm -f $(SOPHYAOBJP)FSAppIrrSmpl.o $(SOPHYAOBJP)cimage.o $(SOPHYAOBJP)cspline.o $(SOPHYAOBJP)datatypes.o $(SOPHYAOBJP)dates.o $(SOPHYAOBJP)difeq.o $(SOPHYAOBJP)dynccd.o $(SOPHYAOBJP)fct1dfit.o $(SOPHYAOBJP)fct2dfit.o $(SOPHYAOBJP)fftmserver.o $(SOPHYAOBJP)fftpserver.o $(SOPHYAOBJP)fftservintf.o $(SOPHYAOBJP)functab.o $(SOPHYAOBJP)generaldata.o $(SOPHYAOBJP)generalfit.o $(SOPHYAOBJP)generalfunc.o $(SOPHYAOBJP)imageop.o $(SOPHYAOBJP)integ.o $(SOPHYAOBJP)median.o $(SOPHYAOBJP)ntoolsinit.o $(SOPHYAOBJP)objfitter.o $(SOPHYAOBJP)perandom.o $(SOPHYAOBJP)poly.o $(SOPHYAOBJP)rk4cdifeq.o $(SOPHYAOBJP)simplesort.o $(SOPHYAOBJP)simplex.o $(SOPHYAOBJP)tabmath.o $(SOPHYAOBJP)toeplitzMatrix.o $(SOPHYAOBJP)datime.o $(SOPHYAOBJP)fftmayer_r4.o $(SOPHYAOBJP)fftmayer_r8.o $(SOPHYAOBJP)fftpackc.o $(SOPHYAOBJP)matxop.o $(SOPHYAOBJP)nbmath.o $(SOPHYAOBJP)nbtri.o
     6         rm -f $(SOPHYAOBJP)FSAppIrrSmpl.o $(SOPHYAOBJP)cimage.o $(SOPHYAOBJP)cspline.o $(SOPHYAOBJP)datatypes.o $(SOPHYAOBJP)dates.o $(SOPHYAOBJP)difeq.o $(SOPHYAOBJP)dynccd.o $(SOPHYAOBJP)fct1dfit.o $(SOPHYAOBJP)fct2dfit.o $(SOPHYAOBJP)fftmserver.o $(SOPHYAOBJP)fftpserver.o $(SOPHYAOBJP)fftservintf.o $(SOPHYAOBJP)functab.o $(SOPHYAOBJP)generaldata.o $(SOPHYAOBJP)generalfit.o $(SOPHYAOBJP)generalfunc.o $(SOPHYAOBJP)imageop.o $(SOPHYAOBJP)integ.o $(SOPHYAOBJP)median.o $(SOPHYAOBJP)ntoolsinit.o $(SOPHYAOBJP)objfitter.o $(SOPHYAOBJP)perandom.o $(SOPHYAOBJP)poly.o $(SOPHYAOBJP)rk4cdifeq.o $(SOPHYAOBJP)simplesort.o $(SOPHYAOBJP)simplex.o $(SOPHYAOBJP)slinparbuff.o $(SOPHYAOBJP)tabmath.o $(SOPHYAOBJP)toeplitzMatrix.o $(SOPHYAOBJP)datime.o $(SOPHYAOBJP)fftmayer_r4.o $(SOPHYAOBJP)fftmayer_r8.o $(SOPHYAOBJP)fftpackc.o $(SOPHYAOBJP)matxop.o $(SOPHYAOBJP)nbmath.o $(SOPHYAOBJP)nbtri.o
    77        rm -f $(SOPHYALIBP)libNTools.a
    88 
    9 $(SOPHYALIBP)libNTools.a : $(SOPHYAOBJP)cimage.o $(SOPHYAOBJP)cspline.o $(SOPHYAOBJP)datatypes.o $(SOPHYAOBJP)dates.o $(SOPHYAOBJP)datime.o $(SOPHYAOBJP)difeq.o $(SOPHYAOBJP)dynccd.o $(SOPHYAOBJP)fct1dfit.o $(SOPHYAOBJP)fct2dfit.o $(SOPHYAOBJP)fftmayer_r4.o $(SOPHYAOBJP)fftmayer_r8.o $(SOPHYAOBJP)fftmserver.o $(SOPHYAOBJP)fftpackc.o $(SOPHYAOBJP)fftpserver.o $(SOPHYAOBJP)fftservintf.o $(SOPHYAOBJP)FSAppIrrSmpl.o $(SOPHYAOBJP)functab.o $(SOPHYAOBJP)generaldata.o $(SOPHYAOBJP)generalfit.o $(SOPHYAOBJP)generalfunc.o $(SOPHYAOBJP)integ.o $(SOPHYAOBJP)matxop.o $(SOPHYAOBJP)median.o $(SOPHYAOBJP)nbmath.o $(SOPHYAOBJP)nbtri.o $(SOPHYAOBJP)ntoolsinit.o $(SOPHYAOBJP)objfitter.o $(SOPHYAOBJP)perandom.o $(SOPHYAOBJP)poly.o $(SOPHYAOBJP)rk4cdifeq.o $(SOPHYAOBJP)simplesort.o $(SOPHYAOBJP)simplex.o $(SOPHYAOBJP)tabmath.o $(SOPHYAOBJP)toeplitzMatrix.o
     9$(SOPHYALIBP)libNTools.a : $(SOPHYAOBJP)cimage.o $(SOPHYAOBJP)cspline.o $(SOPHYAOBJP)datatypes.o $(SOPHYAOBJP)dates.o $(SOPHYAOBJP)datime.o $(SOPHYAOBJP)difeq.o $(SOPHYAOBJP)dynccd.o $(SOPHYAOBJP)fct1dfit.o $(SOPHYAOBJP)fct2dfit.o $(SOPHYAOBJP)fftmayer_r4.o $(SOPHYAOBJP)fftmayer_r8.o $(SOPHYAOBJP)fftmserver.o $(SOPHYAOBJP)fftpackc.o $(SOPHYAOBJP)fftpserver.o $(SOPHYAOBJP)fftservintf.o $(SOPHYAOBJP)FSAppIrrSmpl.o $(SOPHYAOBJP)functab.o $(SOPHYAOBJP)generaldata.o $(SOPHYAOBJP)generalfit.o $(SOPHYAOBJP)generalfunc.o $(SOPHYAOBJP)integ.o $(SOPHYAOBJP)matxop.o $(SOPHYAOBJP)median.o $(SOPHYAOBJP)nbmath.o $(SOPHYAOBJP)nbtri.o $(SOPHYAOBJP)ntoolsinit.o $(SOPHYAOBJP)objfitter.o $(SOPHYAOBJP)perandom.o $(SOPHYAOBJP)poly.o $(SOPHYAOBJP)rk4cdifeq.o $(SOPHYAOBJP)simplesort.o $(SOPHYAOBJP)simplex.o $(SOPHYAOBJP)slinparbuff.o $(SOPHYAOBJP)tabmath.o $(SOPHYAOBJP)toeplitzMatrix.o
    1010        $(AR) $(ARFLAGS) $@ $?
    1111        touch $(SOPHYAINCP)/SophyaConfInfo/libsophya.objlist
     
    3939  $(SOPHYAINCP)matharr.h \
    4040  $(SOPHYAINCP)fioarr.h \
    41   $(SOPHYAINCP)sopemtx.h fftservintf.h \
     41  $(SOPHYAINCP)sopemtx.h \
     42  $(SOPHYAINCP)arrctcast.h fftservintf.h \
    4243  $(SOPHYAINCP)pexceptions.h \
    4344  $(SOPHYAINCP)tmatrix.h \
     
    640641        $(CXXCOMPILE) $(CXXTEMPFLG) -o $@  simplex.cc
    641642 
     643$(SOPHYAOBJP)slinparbuff.o: slinparbuff.cc \
     644  $(SOPHYAINCP)machdefs.h \
     645  $(SOPHYAINCP)pexceptions.h \
     646  $(SOPHYAINCP)machdefs.h slinparbuff.h
     647        $(CXXCOMPILE) $(CXXTEMPFLG) -o $@  slinparbuff.cc
     648 
    642649$(SOPHYAOBJP)tabmath.o: tabmath.cc $(SOPHYAINCP)sopnamsp.h \
    643650  $(SOPHYAINCP)machdefs.h tabmath.h peida.h \
     
    677684  $(SOPHYAINCP)matharr.h \
    678685  $(SOPHYAINCP)fioarr.h \
    679   $(SOPHYAINCP)sopemtx.h fftservintf.h \
     686  $(SOPHYAINCP)sopemtx.h \
     687  $(SOPHYAINCP)arrctcast.h fftservintf.h \
    680688  $(SOPHYAINCP)pexceptions.h \
    681689  $(SOPHYAINCP)tmatrix.h \
  • trunk/SophyaLib/NTools/datatypes.cc

    r2615 r3572  
    99
    1010/* Nouvelle-Fonction */
    11 char * DataName(PBaseDataTypes typ)
     11const char * DataName(PBaseDataTypes typ)
    1212{
    1313switch (typ)
     
    3333
    3434/* Nouvelle-Fonction */
    35 char * DataLongName(PBaseDataTypes typ)
     35const char * DataLongName(PBaseDataTypes typ)
    3636{
    3737switch (typ)
  • trunk/SophyaLib/NTools/datatypes.h

    r244 r3572  
    2323inline PBaseDataTypes DataType(r_8)    {return kr_8;}
    2424
    25 char * DataName(PBaseDataTypes);
    26 char * DataLongName(PBaseDataTypes);
     25const char * DataName(PBaseDataTypes);
     26const char * DataLongName(PBaseDataTypes);
    2727int    DataSize(PBaseDataTypes);
    2828
  • trunk/SophyaLib/NTools/dates.cc

    r2808 r3572  
    3939TimeZone::TimeZone()
    4040{
    41   char* p = getenv("ACQ_TZ");
     41  const char* p = getenv("ACQ_TZ");
    4242#if defined(__DECCXX) || defined(__KCC__) || defined(__aCC__)
    4343  if (!p) p = const_cast<char *>("France");
  • trunk/SophyaLib/NTools/functab.cc

    r2618 r3572  
    5353  Dx = (r_8)(XMax-XMin)/(r_8)(NTab-1);
    5454  Tabul = new r_8[NTab+2];
    55   for(int_4 i=0;i<NTab+2;i++) Tabul[i] = 0.;
     55  for(uint_4 i=0;i<NTab+2;i++) Tabul[i] = 0.;
    5656  PInterpL = PInterpP[0] = PInterpP[1] = NULL;
    5757  if(!Tabul) return;
    58   for(int_4 i=0;i<NTab;i++) Tabul[i+1] = func(XMin + i*Dx);
     58  for(uint_4 i=0;i<NTab;i++) Tabul[i+1] = func(XMin + i*Dx);
    5959  Tabul[0] = Tabul[1];
    6060  Tabul[NTab+1] = Tabul[NTab];
     
    106106    if(PInterpL) delete [] PInterpL;
    107107    PInterpL = new r_8[NTab+2];
    108     for(int_4 i=0;i<NTab+1;i++)
     108    for(uint_4 i=0;i<NTab+1;i++)
    109109      PInterpL[i] = (Tabul[i+1]-Tabul[i])/Dx;
    110110    PInterpL[NTab+1] = 0.;
     
    116116    PInterpP[0] = new r_8[NTab+2];
    117117    PInterpP[1] = new r_8[NTab+2];
    118     for(int_4 i=1;i<=NTab;i++) {
     118    for(uint_4 i=1;i<=NTab;i++) {
    119119      PInterpP[0][i] = (Tabul[i+1]+Tabul[i-1]-2.*Tabul[i])/(2.*Dx*Dx);
    120120      PInterpP[1][i] = (Tabul[i+1]-Tabul[i-1])/(2.*Dx);
  • trunk/SophyaLib/NTools/generalfit.cc

    r3205 r3572  
    279279  Pour ecrire les iterations dans le fichier filename
    280280*/
    281 void GeneralFit::WriteStep(char *filename)
     281void GeneralFit::WriteStep(const char *filename)
    282282{
    283283
  • trunk/SophyaLib/NTools/generalfit.h

    r3083 r3572  
    5353  ~GeneralFit();
    5454
    55   void            WriteStep(char *filename = NULL);
     55  void            WriteStep(const char *filename = NULL);
    5656  void            SetDebug(int level = 0);
    5757  void            SetMaxStep(int n = 100);
  • trunk/SophyaLib/NTools/ntools.h

    r1316 r3572  
    1515#include "poly.h"
    1616#include "datatypes.h"
     17#include "slinparbuff.h"
    1718
    1819#include "cimage.h"
  • trunk/SophyaLib/NTools/objlist.list

    r3083 r3572  
    3131simplesort.o
    3232simplex.o
     33slinparbuff.o
    3334tabmath.o
    3435toeplitzMatrix.o
  • trunk/SophyaLib/NTools/pclassids.h

    r2602 r3572  
    3030   ClassId_HistoErr = 0x0204,
    3131   ClassId_NTuple   = 0x0210,
    32    ClassId_XNTuple  = 0x0211,
    3332
    3433   ClassId_GeneralFitData  = 0x0290,
  • trunk/SophyaLib/NTools/poly.cc

    r3235 r3572  
    710710Poly2& Poly2::operator *= (double a)
    711711{
    712   for (uint_4 i=0; i<NElts(); i++) Element(i) *= a;
     712  for (int_4 i=0; i<NElts(); i++) Element(i) *= a;
    713713  return *this;
    714714}
  • trunk/SophyaLib/NTools/simplex.cc

    r2808 r3572  
    366366  int ilo, ihi, inhi;
    367367  int move = 0;
    368   char* smov[6] = { "None", "Reflection", "ReflecExpand", "ContractHigh", "ContractLow", "ExpandHigh" };
     368  const char* smov[6] = { "None", "Reflection", "ReflecExpand", "ContractHigh", "ContractLow", "ExpandHigh" };
    369369  int movcnt[6] = {0,0,0,0,0,0};
    370370
     
    511511int MinZSimplex::StopReason(string& s)
    512512{
    513   char* sr[5] = { "NoConverg, MaxIterReached", "OK, fm<Tol0", "OK, Df/f<Tol1",
     513  const char* sr[5] = { "NoConverg, MaxIterReached", "OK, fm<Tol0", "OK, Df/f<Tol1",
    514514                  "OK, [Df/f max]Iter<Tol2" "Error - Wrong StopReason" };
    515515  int stop = mStop;
  • trunk/SophyaLib/NTools/slinparbuff.h

    r2322 r3572  
    1 #include "machdefs.h"
    2 #include <iostream>
    3 #include <stdlib.h>
    4 #include <stdio.h>
    5 #include <math.h>
    6 #include <algorithm>
    71#include "pexceptions.h"
    82
     
    1711  friend class SLinParBuffMerger;
    1812
    19   //*****************************************************************
    20   SLinParBuff(uint_4 lenbuf,uint_4 nresynch=0
    21              ,r_8 x0=0.,r_8 y0=0.,bool autoxy0=false)
    22   // ---- Createur
    23   {
    24     if(lenbuf==0) throw RangeCheckError("SLinParBuff: lenbuf==0 !");
    25     mLenBuf = lenbuf;
    26     mNResynch = nresynch;
    27     mX = new r_8[mLenBuf]; mY = new r_8[mLenBuf];
    28     mX0 = x0; mY0 = y0; mAutoXY0 = autoxy0;
    29     Reset();
    30   }
     13  SLinParBuff(uint_4 lenbuf,uint_4 nresynch=0,r_8 x0=0.,r_8 y0=0.,bool autoxy0=false);
     14  SLinParBuff(SLinParBuff& slp);
     15  SLinParBuff(void);
     16  virtual ~SLinParBuff();
    3117
    32   SLinParBuff(SLinParBuff& slp)
    33   // ---- Createur par copie
    34   {
    35     mLenBuf = slp.mLenBuf;
    36     mNResynch = slp.mNResynch;
     18  void GetX0Y0(r_8& x0,r_8& y0);
     19  bool GetAutoX0Y0();
     20  void SetX0Y0(r_8 x0=0.,r_8 y0=0.);
     21  void SetAutoX0Y0(bool autoxy0=false);
    3722
    38     mX = mY = NULL;
    39     if(mLenBuf) {
    40       mX = new r_8[mLenBuf]; mY = new r_8[mLenBuf];
    41       for(uint_4 i=0;i<mLenBuf;i++) {mX[i]=slp.mX[i]; mY[i]=slp.mY[i];}
    42     }
    43 
    44     mX0 = slp.mX0; mY0 = slp.mY0; mAutoXY0 = slp.mAutoXY0;
    45     mNCur = slp.mNCur; mIDeb = slp.mIDeb; mIResynch = slp.mIResynch;
    46     mNResynchEff = slp.mNResynchEff; mNPush = slp.mNPush; mNPop = slp.mNPop;
    47 
    48     mSx = slp.mSx; mSy = slp.mSy;
    49     mSx2 = slp.mSx2; mSy2 = slp.mSy2; mSxy = slp.mSxy;
    50     mSx3 = slp.mSx3; mSx2y = slp.mSx2y;
    51     mSx4 = slp.mSx4;
    52   }
    53 
    54   SLinParBuff(void)
    55   // ---- Createur par defaut
    56   {
    57     mLenBuf = 0; mNResynch = 0;
    58     mX = mY = NULL;
    59     mX0 = mY0 = 0.; mAutoXY0 = false;
    60     Reset();
    61   }
    62 
    63   virtual ~SLinParBuff()
    64   // ---- Destructeur
    65   {
    66     if(mX) delete [] mX; if(mY) delete [] mY;
    67   }
    68 
    69   //*****************************************************************
    70   inline void GetX0Y0(r_8& x0,r_8& y0) {x0 = mX0; y0 = mY0;}
    71   inline bool GetAutoX0Y0() {return mAutoXY0;}
    72 
    73   inline void SetX0Y0(r_8 x0=0.,r_8 y0=0.)
    74   // ---- Stabilite numerique, centrage des x,y en mX0,mY0
    75   {
    76     mX0 = x0; mY0 = y0; mAutoXY0 = false;
    77     ReComputeSum();
    78   }
    79 
    80   inline void SetAutoX0Y0(bool autoxy0=false)
    81   // ---- Stabilite numerique, centrage automatique a <x> et <y>
    82   //      a chaque ReComputeSum
    83   {
    84     mAutoXY0 = autoxy0;
    85     ReComputeSum();
    86   }
    87 
    88   inline void Reset(void)
    89   // ---- Reset des sommes et des compteurs
    90   {
    91     mSx=mSy=mSx2=mSy2=mSxy=mSx3=mSx2y=mSx4=0.;
    92     mNCur = mIDeb = mIResynch = mNResynchEff = mNPush = mNPop = 0;
    93     //if(mX) for(uint_4 i=0;i<mLenBuf;i++) mX[i] = 0.;
    94     //if(mY) for(uint_4 i=0;i<mLenBuf;i++) mY[i] = 0.;
    95   }
    96 
    97   //*****************************************************************
    98   inline uint_4 Pop(void)
    99   // ---- Pour enlever la donnee la plus ancienne
    100   {
    101     if(mNCur==0) return mNCur;
    102     r_8 x=mX[mIDeb]-mX0, y=mY[mIDeb]-mY0; r_8 x2 = x*x;
    103     mSx  -= x;     mSy   -= y;
    104     mSx2 -= x2;    mSy2  -= y*y;   mSxy -= x*y;
    105     mSx3 -= x2*x;  mSx2y -= x2*y;
    106     mSx4 -= x2*x2;
    107     mNCur--; mNPop++;
    108     if(mNCur==0) {Reset(); return mNCur;}
    109     mIDeb++; if(mIDeb==mLenBuf) mIDeb=0;
    110     return mNCur;
    111   }
    112 
    113   inline uint_4 Push(r_8 x,r_8 y)
    114   // ---- Pour ajouter une donnee, la donnee la plus ancienne
    115   //      est enlevee si le buffer est deja plein.
    116   {
    117     uint_4 ifill;
    118     if(mNCur==mLenBuf) {
    119       r_8 X=mX[mIDeb]-mX0, Y=mY[mIDeb]-mY0; r_8 X2 = X*X;
    120       mSx  -= X;     mSy   -=Y;
    121       mSx2 -= X2;    mSy2  -=Y*Y;   mSxy -= X*Y;
    122       mSx3 -= X2*X;  mSx2y -=X*X*Y;
    123       mSx4 -= X2*X2;
    124       ifill = mIDeb;
    125       mIDeb++; if(mIDeb==mLenBuf) mIDeb=0;
    126     } else {
    127       ifill = (mIDeb+mNCur)%mLenBuf;
    128       mNCur++;
    129     }
    130     mX[ifill] = x; mY[ifill] = y;
    131     x -= mX0; y -= mY0; r_8 x2 = x*x;
    132     mSx += x; mSy += y;
    133     mSx2 += x2; mSy2 += y*y; mSxy += x*y;
    134     mSx3 += x2*x; mSx2y += x2*y;
    135     mSx4 += x2*x2;
    136     mNPush++;
    137     // Il faut re-synchroniser pour eviter les derives numeriques
    138     mIResynch++; if(mIResynch == mNResynch) ReComputeSum();
    139     return mNCur;
    140   }
    141 
    142   //*****************************************************************
    143   inline uint_4 ReComputeSum(void)
    144   // ---- Pour re-synchroniser (eviter les derives numeriques).
    145   {
    146     if(mNCur==0) return 0;
    147     // Re-centrage automatique a la moyenne demande:
    148     // Attention, mSx = sum(x-mX0) ==> nouvel mX0 = mSx/N + mX0(ancien)
    149     if(mAutoXY0) {mX0 = mSx/(r_8)mNCur + mX0; mY0 = mSy/(r_8)mNCur + mY0;}
    150     mSx=mSy=mSx2=mSy2=mSxy=mSx3=mSx2y=mSx4=0.;
    151     for(uint_4 i=mIDeb;i<mIDeb+mNCur;i++) {
    152       uint_4 ii = i%mLenBuf;
    153       r_8 x=mX[ii]-mX0, y=mY[ii]-mY0; r_8 x2 = x*x;
    154       mSx  += x;     mSy   += y;
    155       mSx2 += x2;    mSy2  += y*y;   mSxy += x*y;
    156       mSx3 += x2*x;  mSx2y += x2*y;
    157       mSx4 += x2*x2;
    158     }
    159     mIResynch=0;
    160     mNResynchEff++;
    161     return mNCur;
    162   }
    163 
    164   //*****************************************************************
    165   // ---- Retourne le nombre de points
     23  void Reset(void);
     24  uint_4 Pop(void);
     25  uint_4 Push(r_8 x,r_8 y);
    16626  inline uint_4 NPoints(void) {return mNCur;}
    16727
    168   inline r_8 Compute(r_8& mean,bool recomputeXi2=false)
    169   // ---- Calcul <y>, Var(y)
    170   //   recomputeXi2=true : recalcule le xi2 (sigma) avec la courbe et les points
    171   {
    172     mean=0.;
    173     if(mNCur==0) return -1.;
    174     // Moyenne
    175     mean  = mSy/(r_8) mNCur;
    176     // Sigma
    177     r_8 sigma;
    178     if(recomputeXi2) {
    179       sigma=0.;
    180       for(uint_4 i=mIDeb;i<mIDeb+mNCur;i++) {
    181         uint_4 ii = i%mLenBuf;
    182         r_8 s = mean - (mY[ii]-mY0);
    183         sigma += s*s;
    184       }
    185       sigma /= (r_8) mNCur;
    186     } else {
    187       sigma = mSy2/(r_8) mNCur - mean*mean;
    188     }
    189     // gestion du decalage
    190     mean += mY0;
     28  r_8 SumX(void) {return mSx+mNCur*mX0;}
     29  r_8 SumY(void) {return mSy+mNCur*mY0;}
    19130
    192     if(sigma>0.) return sqrt(sigma);
    193     if(sigma<0.) return -sqrt(-sigma);
    194     return sigma;
    195   }
     31  uint_4 ReComputeSum(void);
     32  r_8 Compute(r_8& mean,bool recomputeXi2=false);
     33  r_8 Compute(r_8& a0,r_8 &a1,bool recomputeXi2=false);
     34  r_8 Compute(r_8& a0,r_8 &a1,r_8 &a2,bool recomputeXi2=false);
    19635
    197   inline r_8 Compute(r_8& a0,r_8 &a1,bool recomputeXi2=false)
    198   // ---- Calcul y=a0+a1*x et slin=Var(y-(a0+a1*x))=sqrt(<dy^2>)
    199   //   recomputeXi2=true : recalcule le xi2 avec la courbe et les points
    200   {
    201     a0=a1=0.;
    202     if(mNCur==0) return -1.;
    203     // Fit lineaire
    204     r_8 w = mNCur*mSx2 - mSx*mSx;
    205     if(w==0. || mNCur==1) return -2.;
    206     a0 = (mSx2*mSy - mSx*mSxy)/w;
    207     a1 = (mNCur*mSxy - mSx*mSy )/w;
    208     // Sigma par rapport Fit lineaire
    209     // (On a XI2=mNCur*slin**2 dans notre cas ou les erreurs=1)
    210     r_8 slin;
    211     if(recomputeXi2) {
    212       slin=0.;
    213       for(uint_4 i=mIDeb;i<mIDeb+mNCur;i++) {
    214         uint_4 ii = i%mLenBuf;
    215         r_8 s = a0+a1*(mX[ii]-mX0) - (mY[ii]-mY0);
    216         slin += s*s;
    217       }
    218       slin /= (r_8) mNCur;
    219     } else {
    220       slin = (mSy2 +a0*a0*mNCur +a1*a1*mSx2 -2.*a0*mSy -2.*a1*mSxy +2.*a0*a1*mSx)
    221              / (r_8)mNCur;
    222     }
    223     // gestion du decalage y-y0 = a0 + a1*(x-x0)
    224     a0 = mY0 + a0 - a1*mX0;
    225 
    226     if(slin>0.) return sqrt(slin);
    227     if(slin<0.) return -sqrt(-slin);
    228     return slin;
    229   }
    230 
    231   inline r_8 Compute(r_8& a0,r_8 &a1,r_8 &a2,bool recomputeXi2=false)
    232   // ---- Calcul y=a0+a1*x+a2*x^2 et spar=Var(y-(a0+a1*x+a2*x^2))=sqrt(<dy^2>)
    233   //   recomputeXi2=true : recalcule le xi2 avec la courbe et les points
    234   {
    235     a0=a1=a2=0.;
    236     if(mNCur==0) return -1.;
    237     // Fit parabolique
    238     r_8 w = mSx2*(mSx2*mSx2-mSx3*mSx) -mSx*(mSx3*mSx2-mSx4*mSx) +mNCur*(mSx3*mSx3-mSx4*mSx2);
    239     if(w==0. || mNCur<=2) return -2.;
    240     a2 =  (mSy*(mSx2*mSx2-mSx3*mSx)  - mSxy*(mSx*mSx2-mSx3*mNCur)  + mSx2y*(mSx*mSx-mSx2*mNCur) )/w;
    241     a1 = -(mSy*(mSx3*mSx2-mSx4*mSx)  - mSxy*(mSx2*mSx2-mSx4*mNCur) + mSx2y*(mSx2*mSx-mSx3*mNCur))/w;
    242     a0 =  (mSy*(mSx3*mSx3-mSx4*mSx2) - mSxy*(mSx2*mSx3-mSx4*mSx)   + mSx2y*(mSx2*mSx2-mSx3*mSx) )/w;
    243     // Sigma par rapport Fit parabolique
    244     // (On a XI2=mNCur*spar**2 dans notre cas ou les erreurs=1)
    245     // Le calcul direct du Xi2 n'est pas precis
    246     r_8 spar;
    247     if(recomputeXi2) {
    248       spar=0.;
    249       for(uint_4 i=mIDeb;i<mIDeb+mNCur;i++) {
    250         uint_4 ii = i%mLenBuf;
    251         r_8 s = a0+(a1+a2*(mX[ii]-mX0))*(mX[ii]-mX0) - (mY[ii]-mY0);
    252         spar += s*s;
    253       }
    254       spar /= (r_8) mNCur;
    255     } else {
    256       spar = (mSy2 +a0*a0*mNCur +a1*a1*mSx2 +a2*a2*mSx4 -2.*mSy*a0 -2.*a1*mSxy
    257               -2.*a2*mSx2y +2.*a0*a1*mSx +2.*a0*a2*mSx2 +2.*a1*a2*mSx3)
    258              / (r_8) mNCur;
    259     }
    260     // gestion du decalage y-y0 = a0 + a1*(x-x0) + a2*(x-x0)^2
    261     a0 = mY0 + a0 - a1*mX0 + a2*mX0*mX0;
    262     a1 = a1 - 2.*a2*mX0;
    263 
    264     if(spar>0.) return sqrt(spar);
    265     if(spar<0.) return -sqrt(-spar);
    266     return spar;
    267   }
    268 
    269   //*****************************************************************
    270   void Print(int lp=0)
    271   // ---- Print de l'etat de la classe
    272   //      lp = 0 : parametres
    273   //           1 :            + sommes
    274   //           2 :                     + tableaux
    275   {
    276     cout<<"SLinParBuff(LenBuf="<<mLenBuf<<",NResynch="<<mNResynch
    277         <<",auto_xy0="<<mAutoXY0
    278         <<"): mX0="<<mX0<<" mY0="<<mY0<<" mNCur="<<mNCur<<endl
    279         <<"           NPush="<<mNPush<<" NPop="<<mNPop
    280         <<" NSynchro="<<mNResynchEff<<endl;
    281     if(mNCur==0) return;
    282     if(lp>=2) {
    283       cout<<"X:";
    284       if(mX) for(uint_4 i=0;i<mNCur;i++) {
    285         uint_4 ii = (mIDeb+i)%mLenBuf;
    286         cout<<" "<<mX[ii]; if(i%10==9 || i==mNCur-1) cout<<endl;
    287       }
    288       if(mY) cout<<"Y:";
    289       for(uint_4 i=0;i<mNCur;i++) {
    290         uint_4 ii = (mIDeb+i)%mLenBuf;
    291         cout<<" "<<mY[ii]; if(i%10==9 || i==mNCur-1) cout<<endl;
    292       }
    293     }
    294     if(lp>=1) {
    295       cout<<"...IDeb="<<mIDeb<<" IResynch="<<mIResynch<<endl;
    296       cout<<"   Sx="<<mSx<<" Sx2="<<mSx2<<" Sx3="<<mSx3<<" Sx4="<<mSx4<<endl
    297           <<"   Sy="<<mSy<<" Sy2="<<mSy2<<" Sxy="<<mSxy<<" Sx2y="<<mSx2y<<endl;
    298     }
    299   }
    300 
    301   void PrintCompute(int lp=0)
    302   // ---- Print des valeurs numeriques calculees
    303   //      lp = 0 : valeurs calculees
    304   //           1 :             + sommes
    305   {
    306     bool recompute = (mX==NULL || mY==NULL) ? false: true;
    307     r_8 mean;
    308     r_8 sigma = Compute(mean,recompute), sigmar = Compute(mean,false);
    309     cout<<"SLinParBuff: n="<<NPoints()<<" mean="<<mean
    310         <<" sigma="<<sigma<<" (raw="<<sigmar<<")"<<endl;
    311     r_8 a0,a1,a2;
    312     r_8 slin  = Compute(a0,a1,recompute), slinr  = Compute(a0,a1,false);
    313     cout<<"             a0="<<a0<<" a1="<<a1
    314         <<" slin="<<slin<<" (raw="<<slinr<<")"<<endl;
    315     r_8 spar = Compute(a0,a1,a2,recompute), sparr = Compute(a0,a1,a2,false);
    316     cout<<"             a0="<<a0<<" a1="<<a1<<" a2="<<a2
    317         <<" spar="<<spar<<" (raw="<<sparr<<")"<<endl;
    318     if(lp<1) return;
    319     cout<<"...Sx="<<mSx<<" Sx2="<<mSx2<<" Sx3="<<mSx3<<" Sx4="<<mSx4<<endl
    320         <<"   Sy="<<mSy<<" Sy2="<<mSy2<<" Sxy="<<mSxy<<" Sx2y="<<mSx2y<<endl;
    321   }
     36  void Print(int lp=0);
     37  void PrintCompute(int lp=0);
    32238
    32339protected:
     
    32945  uint_4 mNResynchEff, mNPush, mNPop;
    33046};
     47
    33148
    33249///////////////////////////////////////////////////////////////////
     
    34360class SLinParBuffMerger {
    34461public:
    345   SLinParBuffMerger(void) {mFirst=true; mSlp.Reset();}
    346   SLinParBuffMerger(SLinParBuff& s,bool recompute=false)
    347                    {mFirst=true; mSlp.Reset(); Add(s,recompute);}
    348   virtual ~SLinParBuffMerger(void) {}
    349   inline uint_4  NPoints(void) {return mSlp.NPoints();}
     62  SLinParBuffMerger(void);
     63  SLinParBuffMerger(SLinParBuff& s,bool recompute=false);
     64  virtual ~SLinParBuffMerger(void);
     65
     66  inline uint_4 NPoints(void) {return mSlp.NPoints();}
    35067  inline void Reset(void) {mSlp.Reset(); mFirst=true;}
    351   inline void Add(SLinParBuff& s,bool recompute=false)
    352   {
    353     bool changex0y0=false, AutoXY0_Save; r_8 X0_Save,Y0_Save;
    354     if(mFirst) {
    355       // Cas ou c'est le premier SLinParBuff additionne.
    356       mSlp.mX0=s.mX0; mSlp.mY0=s.mY0; mFirst=false;
    357     } else if(mSlp.mX0!=s.mX0 || mSlp.mY0!=s.mY0) {
    358       // Attention: pour merger il faut avoir les memes mX0,mY0
    359       changex0y0=true;
    360       X0_Save=s.mX0; Y0_Save=s.mY0; AutoXY0_Save=s.mAutoXY0;
    361       s.mX0=mSlp.mX0; s.mY0=mSlp.mY0; s.mAutoXY0=false;
    362       recompute=true;
    363     }
    364     if(recompute) s.ReComputeSum();
    365     mSlp.mNCur += s.mNCur;
    366     mSlp.mSx   += s.mSx;
    367     mSlp.mSy   += s.mSy;
    368     mSlp.mSx2  += s.mSx2;
    369     mSlp.mSy2  += s.mSy2;
    370     mSlp.mSxy  += s.mSxy;
    371     mSlp.mSx3  += s.mSx3;
    372     mSlp.mSx2y += s.mSx2y;
    373     mSlp.mSx4  += s.mSx4;
    374     // Dans le cas ou on a change les X0,Y0, on remet en etat
    375     if(changex0y0) {
    376       s.mX0=X0_Save; s.mY0=Y0_Save;
    377       s.ReComputeSum();
    378       s.mAutoXY0=AutoXY0_Save;
    379     }
    380   }
    381   inline r_8 Compute(r_8& mean)
    382               {return mSlp.Compute(mean);}
    383   inline r_8 Compute(r_8& a0,r_8 &a1)
    384               {return mSlp.Compute(a0,a1);}
    385   inline r_8 Compute(r_8& a0,r_8 &a1,r_8 &a2)
    386               {return mSlp.Compute(a0,a1,a2);}
     68  void Add(SLinParBuff& s,bool recompute=false);
     69  inline r_8 Compute(r_8& mean) {return mSlp.Compute(mean);}
     70  inline r_8 Compute(r_8& a0,r_8 &a1) {return mSlp.Compute(a0,a1);}
     71  inline r_8 Compute(r_8& a0,r_8 &a1,r_8 &a2) {return mSlp.Compute(a0,a1,a2);}
    38772  inline void Print(int lp=0) {mSlp.Print(lp);}
    38873  inline void PrintCompute(int lp=0) {mSlp.PrintCompute(lp);}
  • trunk/SophyaLib/NTools/smakefile

    r3503 r3572  
    44 
    55clean:
    6          rm -f $(SOPHYAOBJP)FSAppIrrSmpl.o $(SOPHYAOBJP)cimage.o $(SOPHYAOBJP)cspline.o $(SOPHYAOBJP)datatypes.o $(SOPHYAOBJP)dates.o $(SOPHYAOBJP)difeq.o $(SOPHYAOBJP)dynccd.o $(SOPHYAOBJP)fct1dfit.o $(SOPHYAOBJP)fct2dfit.o $(SOPHYAOBJP)fftmserver.o $(SOPHYAOBJP)fftpserver.o $(SOPHYAOBJP)fftservintf.o $(SOPHYAOBJP)functab.o $(SOPHYAOBJP)generaldata.o $(SOPHYAOBJP)generalfit.o $(SOPHYAOBJP)generalfunc.o $(SOPHYAOBJP)imageop.o $(SOPHYAOBJP)integ.o $(SOPHYAOBJP)median.o $(SOPHYAOBJP)ntoolsinit.o $(SOPHYAOBJP)objfitter.o $(SOPHYAOBJP)perandom.o $(SOPHYAOBJP)poly.o $(SOPHYAOBJP)rk4cdifeq.o $(SOPHYAOBJP)simplesort.o $(SOPHYAOBJP)simplex.o $(SOPHYAOBJP)tabmath.o $(SOPHYAOBJP)toeplitzMatrix.o $(SOPHYAOBJP)datime.o $(SOPHYAOBJP)fftmayer_r4.o $(SOPHYAOBJP)fftmayer_r8.o $(SOPHYAOBJP)fftpackc.o $(SOPHYAOBJP)matxop.o $(SOPHYAOBJP)nbmath.o $(SOPHYAOBJP)nbtri.o
     6         rm -f $(SOPHYAOBJP)FSAppIrrSmpl.o $(SOPHYAOBJP)cimage.o $(SOPHYAOBJP)cspline.o $(SOPHYAOBJP)datatypes.o $(SOPHYAOBJP)dates.o $(SOPHYAOBJP)difeq.o $(SOPHYAOBJP)dynccd.o $(SOPHYAOBJP)fct1dfit.o $(SOPHYAOBJP)fct2dfit.o $(SOPHYAOBJP)fftmserver.o $(SOPHYAOBJP)fftpserver.o $(SOPHYAOBJP)fftservintf.o $(SOPHYAOBJP)functab.o $(SOPHYAOBJP)generaldata.o $(SOPHYAOBJP)generalfit.o $(SOPHYAOBJP)generalfunc.o $(SOPHYAOBJP)imageop.o $(SOPHYAOBJP)integ.o $(SOPHYAOBJP)median.o $(SOPHYAOBJP)ntoolsinit.o $(SOPHYAOBJP)objfitter.o $(SOPHYAOBJP)perandom.o $(SOPHYAOBJP)poly.o $(SOPHYAOBJP)rk4cdifeq.o $(SOPHYAOBJP)simplesort.o $(SOPHYAOBJP)simplex.o $(SOPHYAOBJP)slinparbuff.o $(SOPHYAOBJP)tabmath.o $(SOPHYAOBJP)toeplitzMatrix.o $(SOPHYAOBJP)datime.o $(SOPHYAOBJP)fftmayer_r4.o $(SOPHYAOBJP)fftmayer_r8.o $(SOPHYAOBJP)fftpackc.o $(SOPHYAOBJP)matxop.o $(SOPHYAOBJP)nbmath.o $(SOPHYAOBJP)nbtri.o
    77        rm -f $(SOPHYALIBP)libNTools.a
    88 
    9 $(SOPHYALIBP)libNTools.a : $(SOPHYAOBJP)cimage.o $(SOPHYAOBJP)cspline.o $(SOPHYAOBJP)datatypes.o $(SOPHYAOBJP)dates.o $(SOPHYAOBJP)datime.o $(SOPHYAOBJP)difeq.o $(SOPHYAOBJP)dynccd.o $(SOPHYAOBJP)fct1dfit.o $(SOPHYAOBJP)fct2dfit.o $(SOPHYAOBJP)fftmayer_r4.o $(SOPHYAOBJP)fftmayer_r8.o $(SOPHYAOBJP)fftmserver.o $(SOPHYAOBJP)fftpackc.o $(SOPHYAOBJP)fftpserver.o $(SOPHYAOBJP)fftservintf.o $(SOPHYAOBJP)FSAppIrrSmpl.o $(SOPHYAOBJP)functab.o $(SOPHYAOBJP)generaldata.o $(SOPHYAOBJP)generalfit.o $(SOPHYAOBJP)generalfunc.o $(SOPHYAOBJP)integ.o $(SOPHYAOBJP)matxop.o $(SOPHYAOBJP)median.o $(SOPHYAOBJP)nbmath.o $(SOPHYAOBJP)nbtri.o $(SOPHYAOBJP)ntoolsinit.o $(SOPHYAOBJP)objfitter.o $(SOPHYAOBJP)perandom.o $(SOPHYAOBJP)poly.o $(SOPHYAOBJP)rk4cdifeq.o $(SOPHYAOBJP)simplesort.o $(SOPHYAOBJP)simplex.o $(SOPHYAOBJP)tabmath.o $(SOPHYAOBJP)toeplitzMatrix.o
     9$(SOPHYALIBP)libNTools.a : $(SOPHYAOBJP)cimage.o $(SOPHYAOBJP)cspline.o $(SOPHYAOBJP)datatypes.o $(SOPHYAOBJP)dates.o $(SOPHYAOBJP)datime.o $(SOPHYAOBJP)difeq.o $(SOPHYAOBJP)dynccd.o $(SOPHYAOBJP)fct1dfit.o $(SOPHYAOBJP)fct2dfit.o $(SOPHYAOBJP)fftmayer_r4.o $(SOPHYAOBJP)fftmayer_r8.o $(SOPHYAOBJP)fftmserver.o $(SOPHYAOBJP)fftpackc.o $(SOPHYAOBJP)fftpserver.o $(SOPHYAOBJP)fftservintf.o $(SOPHYAOBJP)FSAppIrrSmpl.o $(SOPHYAOBJP)functab.o $(SOPHYAOBJP)generaldata.o $(SOPHYAOBJP)generalfit.o $(SOPHYAOBJP)generalfunc.o $(SOPHYAOBJP)integ.o $(SOPHYAOBJP)matxop.o $(SOPHYAOBJP)median.o $(SOPHYAOBJP)nbmath.o $(SOPHYAOBJP)nbtri.o $(SOPHYAOBJP)ntoolsinit.o $(SOPHYAOBJP)objfitter.o $(SOPHYAOBJP)perandom.o $(SOPHYAOBJP)poly.o $(SOPHYAOBJP)rk4cdifeq.o $(SOPHYAOBJP)simplesort.o $(SOPHYAOBJP)simplex.o $(SOPHYAOBJP)slinparbuff.o $(SOPHYAOBJP)tabmath.o $(SOPHYAOBJP)toeplitzMatrix.o
    1010        $(AR) $(ARFLAGS) $@ $?
    1111        touch $(SOPHYAINCP)/SophyaConfInfo/libsophya.objlist
     
    641641        $(CXXCOMPILE) $(CXXTEMPFLG) -o $@  simplex.cc
    642642 
     643$(SOPHYAOBJP)slinparbuff.o: slinparbuff.cc \
     644  $(SOPHYAINCP)machdefs.h \
     645  $(SOPHYAINCP)pexceptions.h \
     646  $(SOPHYAINCP)machdefs.h slinparbuff.h
     647        $(CXXCOMPILE) $(CXXTEMPFLG) -o $@  slinparbuff.cc
     648 
    643649$(SOPHYAOBJP)tabmath.o: tabmath.cc $(SOPHYAINCP)sopnamsp.h \
    644650  $(SOPHYAINCP)machdefs.h tabmath.h peida.h \
Note: See TracChangeset for help on using the changeset viewer.