Index: trunk/SophyaLib/NTools/matrix.cc
===================================================================
--- trunk/SophyaLib/NTools/matrix.cc	(revision 286)
+++ trunk/SophyaLib/NTools/matrix.cc	(revision 288)
@@ -1,3 +1,3 @@
-// $Id: matrix.cc,v 1.2 1999-04-22 16:18:41 ansari Exp $
+// $Id: matrix.cc,v 1.3 1999-05-03 16:55:21 ansari Exp $
 
 #include "machdefs.h"
@@ -10,4 +10,5 @@
 #include "cvector.h"
 #include "generalfit.h"
+#include "tmatrix.h"
 
 // Sous Linux, LN_MAXDOUBLE LN_MINDOUBLE ne semblent pas etre definies
@@ -41,5 +42,5 @@
 //	Construit une matrice 1x1 (pour ppersist).
 //--
-: nr(1), nc(1), ndata(1), nalloc(1), data(new double[1])
+: nr(1), nc(1), ndata(1), nalloc(1), data(new double[1]), FromTmatrix(false)
 {
   // if (r<=0 || c<=0) THROW(rangeCheckErr);
@@ -53,5 +54,5 @@
 //	Construit une matrice de r lignes et c colonnes.
 //--
-: nr(r), nc(c), ndata(r*c), nalloc(r*c), data(new double[r*c])
+: nr(r), nc(c), ndata(r*c), nalloc(r*c), data(new double[r*c]), FromTmatrix(false)
 {
   // if (r<=0 || c<=0) THROW(rangeCheckErr);
@@ -66,5 +67,5 @@
 //	le tableau des valeurs : pas d'allocation.
 //--
-: nr(r), nc(c), ndata(r*c), nalloc(0), data(values)
+: nr(r), nc(c), ndata(r*c), nalloc(0), data(values), FromTmatrix(false)
 {
   END_CONSTRUCTOR
@@ -77,5 +78,5 @@
 //--
 : nr(a.nr), nc(a.nc), ndata(a.nr*a.nc), nalloc(a.nr*a.nc),
-  data(new double[a.nr*a.nc])
+  data(new double[a.nr*a.nc]), FromTmatrix(false)
 {
   memcpy(data, a.data, nalloc * sizeof(double));
@@ -83,9 +84,23 @@
 }
 
+//++
+Matrix::Matrix(const TMatrix<r_8>& a)
+//
+//	Constructeur par copie a partir d'une TMatrix<r_8>.
+//	Attention, les donnees sont partagees.
+//--
+: nr(a.NRows()), nc(a.NCols()), ndata(a.NRows()*a.NCols())
+, nalloc(a.NRows()*a.NCols())
+, data(const_cast<r_8 *>(a.Data()))
+, FromTmatrix(true)
+{
+  END_CONSTRUCTOR
+}
+
 
 Matrix::~Matrix()
 {
   DBASSERT(ndata == nr*nc);
-  if (nalloc) delete[] data;
+  if (nalloc && !FromTmatrix) delete[] data;
 }
 
@@ -102,6 +117,5 @@
 {
   DBASSERT(ndata == nr*nc);
-  for (int i=0; i<ndata; i++)
-    data[i] = 0;
+  for (int i=0; i<ndata; i++) data[i] = 0;
 }
 
Index: trunk/SophyaLib/NTools/matrix.h
===================================================================
--- trunk/SophyaLib/NTools/matrix.h	(revision 286)
+++ trunk/SophyaLib/NTools/matrix.h	(revision 288)
@@ -12,4 +12,6 @@
 
 class GeneralFit;
+#include "tmatrix.h"
+//cmv template <class T> class TMatrix;
 
 // <summary> Operations matricielles </summary>
@@ -18,4 +20,5 @@
   friend class Vector;
   friend class MatrixRC;
+  friend class TMatrix<r_8>;
 public:
   Matrix();
@@ -26,4 +29,6 @@
   // Constructeur par copie
   Matrix(const Matrix& a);
+  // Constructeur par copie a partir d'une TMatrix<r_8>
+  Matrix(const TMatrix<r_8>& a);
   // Destructeur
   virtual ~Matrix();
@@ -151,4 +156,5 @@
   int nalloc;
   r_8* data;
+  bool FromTmatrix;
 };
 
Index: trunk/SophyaLib/NTools/tmatrix.cc
===================================================================
--- trunk/SophyaLib/NTools/tmatrix.cc	(revision 286)
+++ trunk/SophyaLib/NTools/tmatrix.cc	(revision 288)
@@ -1,3 +1,3 @@
-// $Id: tmatrix.cc,v 1.2 1999-04-30 11:02:52 ansari Exp $
+// $Id: tmatrix.cc,v 1.3 1999-05-03 16:55:22 ansari Exp $
 //                         C.Magneville          04/99
 #include "machdefs.h"
@@ -13,5 +13,5 @@
 
 ////////////////////////////////////////////////////////////////
-//************ Createur, Destructeur, gestion des donnees
+//**** Createur, Destructeur
 
 template <class T>
@@ -39,5 +39,5 @@
 template <class T>
 TMatrix<T>::TMatrix(const TMatrix<T>& a)
-// Constructeur par copie (partage des donnees).
+// Constructeur par copie (partage si "a" temporaire).
 : mNr(a.mNr), mNc(a.mNc), mNDBlock(a.mNDBlock)
 {
@@ -46,5 +46,5 @@
 template <class T>
 TMatrix<T>::TMatrix(const TMatrix<T>& a,bool share)
-// Constructeur par copie.
+// Constructeur par copie avec possibilite de forcer le partage ou non.
 : mNr(a.mNr), mNc(a.mNc), mNDBlock(a.mNDBlock,share)
 {
@@ -58,51 +58,16 @@
 
 ////////////////////////////////////////////////////////////////
-template <class T>
-void TMatrix<T>::Clone(const TMatrix<T>& a)
-// Clone (copie de donnees) a partir de "a"
-// Partage des donnees si "a" est temporaire.
-{
-mNDBlock.Clone(a.mNDBlock); mNr = a.mNr; mNc = a.mNc;
-}
-
-template <class T>
-void TMatrix<T>::Reset(T v)
-// Reset de la matrice a "v"
-{
-mNDBlock.Reset(v);
-}
-
-template <class T>
-void TMatrix<T>::Realloc(uint_4 r,uint_4 c,bool force_alloc)
-// Reallocation de place
-{
-if(r==0 || c==0) throw(SzMismatchError("TMatrix::ReSize r ou c==0\n"));
-if(!force_alloc && mNr==r && mNc==c)  return;
-mNr = r; mNc = c; mNDBlock.ReSize(r*c,force_alloc);
-}
-
-////////////////////////////////////////////////////////////////
-//**** Surcharge de TMatrix=TMatrix; TMatrix=<T> b; 
-
-template <class T>
-TMatrix<T>& TMatrix<T>::operator = (T x)
-// Operateur d'affectation depuis scalaire : identite * scalaire.
-{
-if(mNr!=mNc || mNr==0)
-  throw(SzMismatchError("TMatrix::operator= mNc!=mNr ou ==0\n"));
-for(uint_4 r=0;r<mNr;r++) for(uint_4 c=0;c<mNc;c++)
-                               (*this)(r,c) = (r==c)? x: 0;
-return *this;
-}
-
-template <class T>
-TMatrix<T>& TMatrix<T>::operator = (const TMatrix<T>& a)
-// Operateur d'affectation A=a
-// (Re)allocation de la place memoire
-// sauf si "a" est cree temporairement (partage)
-{
-if(this == &a) return *this;
-Clone(a);
-return *this;
+// Operations matricielles
+template <class T>
+TMatrix<T> TMatrix<T>::Transpose(void) const
+// Transposition
+{
+TMatrix<T> a; a.Clone(*this); a.SetTemp(true);
+a.mNr = mNc; a.mNc = mNr;
+{for(int i=0; i<mNr; i++)
+  for(int j=0; j<mNc; j++) {
+    a(j,i) = (*this)(i,j);
+}}
+return a;
 }
 
@@ -115,5 +80,5 @@
 // Impression de la sous-matrice (i:i+ni-1,i:j+nj-1)
 {
-  os<<"TMatrix::Print("<<mNr<<","<<mNc<<")"<<endl;
+os<<"TMatrix::Print("<<mNr<<","<<mNc<<")"<<endl;
 if(lp>0)
   {os<<"       this="<<this<<endl; mNDBlock.Print(0,0);}
@@ -129,54 +94,5 @@
 
 ////////////////////////////////////////////////////////////////
-//**** Surcharge de +=,-=,*=,/= (INPLACE): TMatrix += <T> b;
-
-template <class T>
-TMatrix<T>& TMatrix<T>::operator += (T b)
-{
-mNDBlock += b;
-return *this;
-}
-
-template <class T>
-TMatrix<T>& TMatrix<T>::operator -= (T b)
-{
-mNDBlock -= b;
-return *this;
-}
-
-template <class T>
-TMatrix<T>& TMatrix<T>::operator *= (T b)
-{
-mNDBlock *= b;
-return *this;
-}
-
-template <class T>
-TMatrix<T>& TMatrix<T>::operator /= (T b)
-{
-mNDBlock /= b;
-return *this;
-}
-
-////////////////////////////////////////////////////////////////
-//**** Surcharge de +=,-=,*= (INPLACE): TMatrix += TMatrix;
-
-template <class T>
-TMatrix<T>& TMatrix<T>::operator += (const TMatrix<T>& a)
-{
-if(mNr==0 || mNc==0 || mNr!=a.mNr || mNc!=a.mNc)
-  throw(SzMismatchError("TMatrix::operator+=A size mismatch"));
-mNDBlock += a.mNDBlock;
-return *this;
-}
-
-template <class T>
-TMatrix<T>& TMatrix<T>::operator -= (const TMatrix<T>& a)
-{
-if(mNr==0 || mNc==0 || mNr!=a.mNr || mNc!=a.mNc)
-  throw(SzMismatchError("TMatrix::operator-=A size mismatch"));
-mNDBlock -= a.mNDBlock;
-return *this;
-}
+//**** Surcharge de *= (INPLACE): TMatrix *= TMatrix;
 
 template <class T>
@@ -209,47 +125,30 @@
 
 ////////////////////////////////////////////////////////////////
-//**** Surcharge de +,-,*: TMatrix = TMatrix + TMatrix;
-
-template <class T>
-TMatrix<T> TMatrix<T>::Add(const TMatrix<T>& b) const
-// C = A(this)+B
+//**** Pour surcharge d'operateurs C = A (+,-,*,/) B
+
+template <class T> TMatrix<T> TMatrix<T>::Add(const TMatrix<T>& b) const
 {
 if(mNr!=b.mNr || mNc!=b.mNc)
-  throw(SzMismatchError("NDataBlock operator C=A+B size mismatch\n"));
-TMatrix<T> result; result.SetTemp(true);
-if(b.IsTemp()) {
-  result.Clone(b);
-  result += *this;
-} else {
-  result.Clone(*this);
-  result += b;
-}
+    throw(SzMismatchError("NDataBlock operator C=A+B size mismatch\n"));
+TMatrix<T> result; result.SetTemp(true); result.mNr=mNr; result.mNc=mNc;
+result.mNDBlock = mNDBlock+b.mNDBlock;
 return result;
 }
 
-template <class T>
-TMatrix<T> TMatrix<T>::Sub(const TMatrix<T>& b) const
-// C = A(this)-B
+template <class T> TMatrix<T> TMatrix<T>::Sub(const TMatrix<T>& b) const
 {
 if(mNr!=b.mNr || mNc!=b.mNc)
-  throw(SzMismatchError("NDataBlock operator C=A-B size mismatch\n"));
-TMatrix<T> result; result.SetTemp(true);
-if(b.IsTemp()) {
-  result.Clone(b);
-  result.mNDBlock = (*this).mNDBlock - result.mNDBlock;
-} else {
-  result.Clone(*this);
-  result -= b;
-}
+    throw(SzMismatchError("NDataBlock operator C=A-B size mismatch\n"));
+TMatrix<T> result; result.SetTemp(true); result.mNr=mNr; result.mNc=mNc;
+result.mNDBlock = mNDBlock-b.mNDBlock;
 return result;
 }
 
-template <class T>
-TMatrix<T> TMatrix<T>::Mul(const TMatrix<T>& b) const
+template <class T> TMatrix<T> TMatrix<T>::Mul(const TMatrix<T>& b) const
 // C = A(this)*B : Cij = Aik Bkj (allocation forcee dans tous les cas)
 {
 if(mNr==0 || mNc==0 || b.mNc==0 || mNc!=b.mNr)
   throw(SzMismatchError("NDataBlock operator C=A*B size mismatch\n"));
-TMatrix<T> r; r.SetTemp(true); r.Realloc(mNr,b.mNc,true);
+TMatrix<T> r; r.SetTemp(true); r.ReSize(mNr,b.mNc);
 T *ai,*aik,*bj,*bkj,*ri,*rij;
 for(ri=const_cast<T *>(r.Data()),ai=const_cast<T *>(Data());
@@ -317,5 +216,5 @@
 is.Get(itab,2);
 if (dobj == NULL) dobj = new TMatrix<T>(itab[0],itab[1]);
-else dobj->Realloc(itab[0],itab[1],false);
+else dobj->ReSize(itab[0],itab[1]);
 // On lit le NDataBlock
 //cmv   encore des problemes avec les consteries
@@ -331,5 +230,5 @@
 uint_4 itab[2];
 itab[0] = dobj->NRows();
-itab[1] = dobj->NCol();
+itab[1] = dobj->NCols();
 os.Put(itab,2);
 // On ecrit le NDataBlock
Index: trunk/SophyaLib/NTools/tmatrix.h
===================================================================
--- trunk/SophyaLib/NTools/tmatrix.h	(revision 286)
+++ trunk/SophyaLib/NTools/tmatrix.h	(revision 288)
@@ -16,4 +16,5 @@
 class TMatrix : public AnyDataObj {
 public:
+
   // Creation / destruction 
   TMatrix();
@@ -26,19 +27,25 @@
   // Temporaire?
   inline bool IsTemp(void) const {return mNDBlock.IsTemp();}
-  inline void SetTemp(bool temp=false) const
-                     {mNDBlock.SetTemp(temp);}
+  inline void SetTemp(bool temp=false) const {mNDBlock.SetTemp(temp);}
 
   // Gestion taille/Remplissage
-  void Clone(const TMatrix<T>& a);
-  void Reset(T v=0);
-  void Realloc(uint_4 r,uint_4 c,bool force_alloc=true);
+  inline void Clone(const TMatrix<T>& a)  // Clone: copie des donnees de "a"
+              {mNDBlock.Clone(a.mNDBlock); mNr = a.mNr; mNc = a.mNc;}
+  inline void Reset(T v=0) {mNDBlock.Reset(v);}
+  inline void ReSize(uint_4 r,uint_4 c)  // Reallocation de place
+       {if(r==0||c==0) throw(SzMismatchError("TMatrix::ReSize r ou c==0\n"));
+        mNr = r; mNc = c; mNDBlock.ReSize(r*c);}
 
   // Informations pointeur/data
   inline int NRows() const {return mNr;}
-  inline int NCol()  const {return mNc;}
-  T const& operator()(uint_4 r,uint_4 c) const
-                     {return *(mNDBlock.Begin()+r*mNc+c);}
-  T&       operator()(uint_4 r,uint_4 c)
-                     {return *(mNDBlock.Begin()+r*mNc+c);}
+  inline int NCols() const {return mNc;}
+  inline T const& operator()(uint_4 r,uint_4 c) const
+                            {return *(mNDBlock.Begin()+r*mNc+c);}
+  inline T&       operator()(uint_4 r,uint_4 c)
+                            {return *(mNDBlock.Begin()+r*mNc+c);}
+  inline T const& operator[](uint_4 ip) const
+                            {return *(mNDBlock.Begin()+ip);}
+  inline T&       operator[](uint_4 ip)
+                            {return *(mNDBlock.Begin()+ip);}
   inline       T* Data()       {return mNDBlock.Begin();}
   inline const T* Data() const {return mNDBlock.Begin();}
@@ -46,7 +53,16 @@
   inline const NDataBlock<T>& DataBlock() const {return mNDBlock;}
 
+  // Operations matricielles
+  TMatrix<T> Transpose(void) const;
+
   // Operateur d'affectation
-  TMatrix<T>& operator = (T x);
-  TMatrix<T>& operator = (const TMatrix<T>& a);
+  // A = x (matrice diagonale x*Identite)
+  inline TMatrix<T>& operator = (T x)
+    {if(mNr!=mNc || mNr==0) throw(SzMismatchError("TMatrix::operator= mNc!=mNr ou ==0\n"));
+     for(uint_4 r=0;r<mNr;r++) for(uint_4 c=0;c<mNc;c++) (*this)(r,c)=(r==c)?x:0;
+     return *this;}
+  // A = B : partage les donnees si "a" est temporaire, clone sinon.
+  inline TMatrix<T>& operator = (const TMatrix<T>& a)
+                {if(this == &a) return *this; CloneOrShare(a); return *this;}
 
   // Impression
@@ -55,15 +71,22 @@
               {Print(cout,lp,i0,ni,j0,nj);}
 
-  // Surcharge d'operateurs INPLACE
-  TMatrix<T>& operator += (T b);
-  TMatrix<T>& operator -= (T b);
-  TMatrix<T>& operator *= (T b);
-  TMatrix<T>& operator /= (T b);
+  // Surcharge d'operateurs INPLACE: A (+=,-=,*=,/=) (T) x
+  inline TMatrix<T>& operator += (T b) {mNDBlock += b; return *this;}
+  inline TMatrix<T>& operator -= (T b) {mNDBlock -= b; return *this;}
+  inline TMatrix<T>& operator *= (T b) {mNDBlock *= b; return *this;}
+  inline TMatrix<T>& operator /= (T b) {mNDBlock /= b; return *this;}
 
-  TMatrix<T>& operator += (const TMatrix<T>& a);
-  TMatrix<T>& operator -= (const TMatrix<T>& a);
+  // Surcharge d'operateurs INPLACE: A (+=,-=,*=,/=) B
+  inline TMatrix<T>& operator += (const TMatrix<T>& a)
+              {if(mNr==0 || mNc==0 || mNr!=a.mNr || mNc!=a.mNc)
+                  throw(SzMismatchError("TMatrix::operator+=A size mismatch"));
+               mNDBlock += a.mNDBlock; return *this;}
+  inline TMatrix<T>& operator -= (const TMatrix<T>& a)
+              {if(mNr==0 || mNc==0 || mNr!=a.mNr || mNc!=a.mNc)
+                  throw(SzMismatchError("TMatrix::operator-=A size mismatch"));
+               mNDBlock -= a.mNDBlock; return *this;}
   TMatrix<T>& operator *= (const TMatrix<T>& a);
 
-  // Pour surcharge d'operateurs
+  // Pour surcharge d'operateurs C = A (+,-,*) B
   TMatrix<T> Add(const TMatrix<T>& b) const;
   TMatrix<T> Sub(const TMatrix<T>& b) const;
@@ -71,75 +94,62 @@
 
 protected:
+  // partage les donnees si "a" temporaire, clone sinon.
+  inline void CloneOrShare(const TMatrix<T>& a)
+              {mNDBlock.CloneOrShare(a.mNDBlock); mNr=a.mNr; mNc=a.mNc;}
+  // Share: partage les donnees de "a"
+  inline void Share(const TMatrix<T>& a)
+              {mNDBlock.Share(a.mNDBlock); mNr=a.mNr; mNc=a.mNc;}
+
   uint_4 mNr,mNc;
   NDataBlock<T> mNDBlock;
 };
 
-/////////////////////////////////////////////////////////////////////////
-template<class T>
+////////////////////////////////////////////////////////////////
+// Impression
+
+template <class T>
 inline ostream& operator << (ostream& os, const TMatrix<T>& a)
-                      {a.Print(os); return(os);}
+                            {a.Print(os); return(os);}
 
-/////////////////////////////////////////////////////////////////////////
-template<class T>
-inline TMatrix<T> operator + (const TMatrix<T>& a, double b) {
-  TMatrix<T> result(a,false); result.DataBlock().SetTemp(true);
-  result += b;
-  return result;
-}
+////////////////////////////////////////////////////////////////
+// Surcharge d'operateurs A (+=,-=,*=,/=) (T) x
 
-template<class T>
-inline TMatrix<T> operator + (double b,const TMatrix<T>& a) {
-  TMatrix<T> result(a,false); result.DataBlock().SetTemp(true);
-  result += b;
-  return result;
-}
+template <class T> inline TMatrix<T> operator + (const TMatrix<T>& a, T b)
+    {TMatrix<T> result(a); result.SetTemp(true); result += b; return result;}
 
-template<class T>
-inline TMatrix<T> operator - (const TMatrix<T>& a, double b) {
-  TMatrix<T> result(a,false); result.DataBlock().SetTemp(true);
-  result -= b;
-  return result;
-}
+template <class T> inline TMatrix<T> operator + (T b,const TMatrix<T>& a)
+    {TMatrix<T> result(a); result.SetTemp(true); result += b; return result;}
 
-template<class T>
-inline TMatrix<T> operator - (double b,const TMatrix<T>& a) {
-  TMatrix<T> result(a,false); result.DataBlock().SetTemp(true);
-  result.DataBlock() = b - result.DataBlock();
-  return result;
-}
+template <class T> inline TMatrix<T> operator - (const TMatrix<T>& a, T b)
+    {TMatrix<T> result(a); result.SetTemp(true); result -= b; return result;}
 
-template<class T>
-inline TMatrix<T> operator * (const TMatrix<T>& a, double b) {
-  TMatrix<T> result(a,false); result.DataBlock().SetTemp(true);
-  result *= b;
-  return result;
-}
+template <class T> inline TMatrix<T> operator - (T b,const TMatrix<T>& a)
+    {TMatrix<T> result(a); result.SetTemp(true);
+     result.DataBlock() = b-result.DataBlock(); return result;}
 
-template<class T>
-inline TMatrix<T> operator * (double b,const TMatrix<T>& a) {
-  TMatrix<T> result(a,false); result.DataBlock().SetTemp(true);
-  result *= b;
-  return result;
-}
+template <class T> inline TMatrix<T> operator * (const TMatrix<T>& a, T b)
+    {TMatrix<T> result(a); result.SetTemp(true); result *= b; return result;}
 
-template<class T>
-inline TMatrix<T> operator / (const TMatrix<T>& a, double b) {
-  TMatrix<T> result(a,false); result.DataBlock().SetTemp(true);
-  result /= b;
-  return result;
-}
+template <class T> inline TMatrix<T> operator * (T b,const TMatrix<T>& a)
+    {TMatrix<T> result(a); result.SetTemp(true); result *= b; return result;}
 
-/////////////////////////////////////////////////////////////////////////
-template<class T>
+template <class T> inline TMatrix<T> operator / (const TMatrix<T>& a, T b)
+    {TMatrix<T> result(a); result.SetTemp(true); result /= b; return result;}
+
+////////////////////////////////////////////////////////////////
+// Surcharge d'operateurs C = A (+,-,*,/) B
+
+template <class T>
 inline TMatrix<T> operator + (const TMatrix<T>& a,const TMatrix<T>& b)
-                             {return a.Add(b);}
+                  {return a.Add(b);}
 
-template<class T>
+template <class T>
 inline TMatrix<T> operator - (const TMatrix<T>& a,const TMatrix<T>& b)
-                             {return a.Sub(b);}
+                  {return a.Sub(b);}
 
-template<class T>
+template <class T>
 inline TMatrix<T> operator * (const TMatrix<T>& a,const TMatrix<T>& b)
-                             {return a.Mul(b);}
+                  {return a.Mul(b);}
+
 
 /////////////////////////////////////////////////////////////////////////
@@ -147,18 +157,15 @@
 template <class T>
 class FIO_TMatrix : public  PPersist  {
-
 public:
-	    FIO_TMatrix();
-	    FIO_TMatrix(string const & filename); 
-	    FIO_TMatrix(const TMatrix<T> & obj);
-	    FIO_TMatrix(TMatrix<T> * obj);
-  virtual   ~FIO_TMatrix();
-
-  virtual   AnyDataObj* DataObj();
+  FIO_TMatrix();
+  FIO_TMatrix(string const & filename); 
+  FIO_TMatrix(const TMatrix<T> & obj);
+  FIO_TMatrix(TMatrix<T> * obj);
+  virtual ~FIO_TMatrix();
+  virtual AnyDataObj* DataObj();
   inline operator TMatrix<T>() { return(*dobj); }
-
 protected :
-  virtual void       ReadSelf(PInPersist&);           
-  virtual void       WriteSelf(POutPersist&) const;  
+  virtual void ReadSelf(PInPersist&);           
+  virtual void WriteSelf(POutPersist&) const;  
   TMatrix<T> * dobj;
   bool ownobj;
