Index: trunk/SophyaLib/BaseTools/ndatablock.cc
===================================================================
--- trunk/SophyaLib/BaseTools/ndatablock.cc	(revision 267)
+++ trunk/SophyaLib/BaseTools/ndatablock.cc	(revision 268)
@@ -1,2 +1,5 @@
+// Gestion de block de donnees avec partage de references
+//                         C.Magneville          04/99
+// LAL (Orsay) / IN2P3-CNRS  DAPNIA/SPP (Saclay) / CEA
 #include "machdefs.h"
 #include <stdio.h>
@@ -7,5 +10,7 @@
 #include "ndatablock.h"
 
-#define DEBUG_NDATABLOCK
+using namespace PlanckDPC;
+
+//define DEBUG_NDATABLOCK
 
 #ifdef DEBUG_NDATABLOCK
@@ -25,4 +30,8 @@
 : mSz(0), mSRef(NULL), mIsTemp(false)
 {
+#ifdef DEBUG_NDATABLOCK
+cout<<"DEBUG_NDataBlock::NDataBlock("<<this<<",n="<<n<<")"<<endl;
+#endif
+
 Alloc(n);
 }
@@ -33,4 +42,9 @@
 : mSz(0), mSRef(NULL), mIsTemp(false)
 {
+#ifdef DEBUG_NDATABLOCK
+cout<<"DEBUG_NDataBlock::NDataBlock("<<this
+    <<",data="<<data<<",br="<<br<<")"<<endl;
+#endif
+
 Alloc(n,data,br);
 }
@@ -41,8 +55,11 @@
 : mSz(0), mSRef(NULL), mIsTemp(false)
 {
-}
-
-template <class T>
-NDataBlock<T>::NDataBlock(NDataBlock<T>& a)
+#ifdef DEBUG_NDATABLOCK
+cout<<"DEBUG_NDataBlock::NDataBlock("<<this<<") default"<<endl;
+#endif
+}
+
+template <class T>
+NDataBlock<T>::NDataBlock(const NDataBlock<T>& a)
 // Createur par copie
 // ATTENTION: partage les donnees avec "a"
@@ -51,12 +68,22 @@
 : mSz(0), mSRef(NULL), mIsTemp(false)
 {
+#ifdef DEBUG_NDATABLOCK
+cout<<"DEBUG_NDataBlock::NDataBlock("<<this<<",&a="<<&a<<") a.(mSz="
+    <<a.mSz<<" mSRef="<<a.mSRef<<" IsTemp="<<a.IsTemp()<<")"<<endl;
+#endif
+
 Share(a);
 }
 
 template <class T>
-NDataBlock<T>::NDataBlock(NDataBlock<T>& a,bool share)
+NDataBlock<T>::NDataBlock(const NDataBlock<T>& a,bool share)
 // Createur avec choix de partager ou non
 : mSz(0), mSRef(NULL), mIsTemp(false)
 {
+#ifdef DEBUG_NDATABLOCK
+cout<<"DEBUG_NDataBlock::NDataBlock("<<this<<",&a="<<&a<<",sh=<<"<<share<<")"
+    <<" a.(mSz="<<a.mSz<<" mSRef="<<a.mSRef<<" IsTemp="<<a.IsTemp()<<")"<<endl;
+#endif
+
 if(share) Share(a); else Clone(a);
 }
@@ -66,16 +93,21 @@
 // Destructeur
 {
+#ifdef DEBUG_NDATABLOCK
+cout<<"DEBUG_NDataBlock::~NDataBlock("<<this<<")"<<endl;
+#endif
+
 Delete();
 }
 
 template <class T>
-void NDataBlock<T>::SetTemp(bool temp)
+void NDataBlock<T>::SetTemp(bool temp) const 
 // Set temporary
 {
 mIsTemp=temp;
-    #ifdef DEBUG_NDATABLOCK
-    cout<<"DEBUG_NDataBlock::SetTemp("<<this<<","<<temp
-        <<"), mSz="<<mSz<<" mSRef="<<mSRef<<" IsTemp="<<mIsTemp<<endl;
-    #endif
+
+#ifdef DEBUG_NDATABLOCK
+cout<<"DEBUG_NDataBlock::SetTemp("<<this<<","<<temp
+    <<"), mSz="<<mSz<<" mSRef="<<mSRef<<" IsTemp="<<mIsTemp<<endl;
+#endif
 }
 
@@ -99,9 +131,10 @@
 //   NDataBlock A(5); A.FillFrom(5,x);
 {
-    #ifdef DEBUG_NDATABLOCK
-    cout<<"DEBUG_NDataBlock::Alloc("<<this<<","
-        <<n<<","<<data<<","<<br<<") mSz="<<mSz
-        <<" mSRef="<<mSRef<<" IsTemp="<<mIsTemp<<endl;
-    #endif
+#ifdef DEBUG_NDATABLOCK
+cout<<"DEBUG_NDataBlock::Alloc("<<this<<","
+    <<n<<","<<data<<","<<br<<") mSz="<<mSz
+    <<" mSRef="<<mSRef<<" IsTemp="<<mIsTemp<<endl;
+#endif
+
 if(n==0) throw(SzMismatchError("NDataBlock::Alloc n==0\n"));
 if(mSRef) Delete();
@@ -111,15 +144,16 @@
 if(data) mSRef->data = data; else mSRef->data = new T[n];
 mSRef->bridge = br;
-    #ifdef DEBUG_NDATABLOCK
-    if(!data) NallocData++; NallocSRef++;
-    cout<<"...DEBUG_NDataBlock::Alloc mSz="<<mSz<<" mSRef="<<mSRef
-        <<" mSRef->nref="<<mSRef->nref<<" mSRef->data="<< mSRef->data
-        <<" mSRef->bridge="<<mSRef->bridge<<" IsTemp="<<mIsTemp
-        <<" Total("<<NallocData<<","<<NallocSRef<<")"<<endl;
-    #endif
-}
-
-template <class T>
-void NDataBlock<T>::Clone(NDataBlock<T>& a)
+
+#ifdef DEBUG_NDATABLOCK
+if(!data) NallocData++; NallocSRef++;
+cout<<"...DEBUG_NDataBlock::Alloc mSz="<<mSz<<" mSRef="<<mSRef
+    <<" mSRef->nref="<<mSRef->nref<<" mSRef->data="<< mSRef->data
+    <<" mSRef->bridge="<<mSRef->bridge<<" IsTemp="<<mIsTemp
+    <<" Total("<<NallocData<<","<<NallocSRef<<")"<<endl;
+#endif
+}
+
+template <class T>
+void NDataBlock<T>::Clone(const NDataBlock<T>& a)
 // Clone (copie de donnee) a partir de "a"
 {
@@ -135,13 +169,14 @@
 
 template <class T>
-void NDataBlock<T>::Share(NDataBlock<T>& a)
+void NDataBlock<T>::Share(const NDataBlock<T>& a)
 // Partage des donnees avec "a"
 {
-    #ifdef DEBUG_NDATABLOCK
-    cout<<"DEBUG_NDataBlock::Share("<<this<<","<<&a<<")";
-    if(&a!=NULL) cout<<" a.(mSz="<<a.mSz<<" mSRef="<<a.mSRef
-                     <<" IsTemp="<<a.IsTemp()<<")";
-    cout<<", mSz="<<mSz<<" mSRef="<<mSRef<<" IsTemp="<<mIsTemp<<endl;
-    #endif
+#ifdef DEBUG_NDATABLOCK
+cout<<"DEBUG_NDataBlock::Share("<<this<<","<<&a<<")";
+if(&a!=NULL) cout<<" a.(mSz="<<a.mSz<<" mSRef="<<a.mSRef
+                 <<" IsTemp="<<a.IsTemp()<<")";
+cout<<", mSz="<<mSz<<" mSRef="<<mSRef<<" IsTemp="<<mIsTemp<<endl;
+#endif
+
 if(&a==NULL) throw(NullPtrError("NDataBlock::Share  &a==NULL\n"));
 // on ne peut partager si "a" pas alloue
@@ -149,9 +184,10 @@
 if(mSRef) Delete();
 mSz = a.mSz; mSRef = a.mSRef; mSRef->nref++;
-    #ifdef DEBUG_NDATABLOCK
-    cout<<"...DEBUG_NDataBlock::Share mSz="<<mSz<<" mSRef="<<mSRef
-        <<" mSRef->nref="<<mSRef->nref<<" mSRef->data="<< mSRef->data
-        <<" mSRef->bridge="<<mSRef->bridge<<" IsTemp="<<mIsTemp<<endl;
-    #endif
+
+#ifdef DEBUG_NDATABLOCK
+cout<<"...DEBUG_NDataBlock::Share mSz="<<mSz<<" mSRef="<<mSRef
+    <<" mSRef->nref="<<mSRef->nref<<" mSRef->data="<< mSRef->data
+    <<" mSRef->bridge="<<mSRef->bridge<<" IsTemp="<<mIsTemp<<endl;
+#endif
 }
 
@@ -160,26 +196,31 @@
 // Pour detruire les pointeurs en tenant compte des references
 {
-    #ifdef DEBUG_NDATABLOCK
-    cout<<"DEBUG_NDataBlock::Delete("<<this<<") mSz="<<mSz
-        <<" mSRef="<<mSRef<<" IsTemp="<<mIsTemp;
-    if(mSRef)
-      cout<<" mSRef->nref="<<mSRef->nref<<" mSRef->data="
-          <<mSRef->data<<" mSRef->bridge="<<mSRef->bridge;
-    cout<<endl;
-    #endif
+#ifdef DEBUG_NDATABLOCK
+cout<<"DEBUG_NDataBlock::Delete("<<this<<") mSz="<<mSz
+    <<" mSRef="<<mSRef<<" IsTemp="<<mIsTemp;
+if(mSRef)
+  cout<<" mSRef->nref="<<mSRef->nref<<" mSRef->data="
+      <<mSRef->data<<" mSRef->bridge="<<mSRef->bridge;
+cout<<endl;
+#endif
+
 if(mSRef==NULL) return; // cas du createur par defaut
 mSRef->nref--;
 if(mSRef->nref != 0) {
-    #ifdef DEBUG_NDATABLOCK
-      cout<<"...DEBUG_NDataBlock::Delete() pas de desallocation il reste nref="
-          <<mSRef->nref<<endl;
-    #endif
+
+#ifdef DEBUG_NDATABLOCK
+cout<<"...DEBUG_NDataBlock::Delete() pas de desallocation il reste nref="
+    <<mSRef->nref<<endl;
+#endif
+
   return;
 }
-    #ifdef DEBUG_NDATABLOCK
-    if(!mSRef->bridge) NallocData--; NallocSRef--;
-    cout<<"...DEBUG_NDataBlock::Delete() desallocation complete il reste nref="
-        <<mSRef->nref<<" Total("<<NallocData<<","<<NallocSRef<<")"<<endl;
-    #endif
+
+#ifdef DEBUG_NDATABLOCK
+if(!mSRef->bridge) NallocData--; NallocSRef--;
+cout<<"...DEBUG_NDataBlock::Delete() desallocation complete il reste nref="
+    <<mSRef->nref<<" Total("<<NallocData<<","<<NallocSRef<<")"<<endl;
+#endif
+
 // Si il y a un Bridge les donnees ne n'appartiennent pas, on detruit le Bridge
 if(mSRef->bridge) delete mSRef->bridge; 
@@ -223,5 +264,5 @@
 
 template <class T>
-void NDataBlock<T>::Print(size_t i1,size_t n)
+void NDataBlock<T>::Print(ostream& os,size_t i1,size_t n) const
 // Impression de n elements a partir de i1
 {
@@ -229,7 +270,7 @@
 T* p = NULL; Bridge* br = NULL;
 if(mSRef) {nr = mSRef->nref; p = mSRef->data; br = mSRef->bridge;}
-cout<<"NDataBlock::Print("<<this<<",Sz="<<mSz<<",IsTemp="<<mIsTemp<<")\n"
-    <<"            mSRef="<<mSRef<<"(nref="<<nr<<",data="<<p
-    <<",bridge="<<br<<")"<<endl;
+os<<"NDataBlock::Print("<<this<<",Sz="<<mSz<<",IsTemp="<<mIsTemp<<")\n"
+  <<"            mSRef="<<mSRef<<"(nref="<<nr<<",data="<<p
+  <<",bridge="<<br<<")"<<endl;
 if(i1>=mSz || n<=0 || !p) return;
 size_t i2 = i1+n; if(i2>mSz) i2=mSz;
@@ -237,8 +278,8 @@
 while(i1<i2) {
   enl = false;
-  cout<<" "<<(*this)(i1);  i1++;
-  if(im==10) {cout<<"\n"; im=1; enl=true;} else im++;
-}
-if(!enl) cout<<endl;
+  os<<" "<<(*this)(i1);  i1++;
+  if(im==8) {os<<"\n"; im=1; enl=true;} else im++;
+}
+if(!enl) os<<endl;
 }
 
@@ -246,9 +287,15 @@
 
 template <class T>
-NDataBlock<T>& NDataBlock<T>::operator = (NDataBlock<T>& a)
+NDataBlock<T>& NDataBlock<T>::operator = (const NDataBlock<T>& a)
 // surcharge avec partage des donnees
 // Ecriture:  NDataBlock a; a = b;
 //            NDataBlock a(10); a = b; (a est re-affecte)
 {
+#ifdef DEBUG_NDATABLOCK
+cout<<"DEBUG_NDataBlock::operator=("<<this<<","<<&a<<") a.(mSz="
+    <<a.mSz<<" mSRef="<<a.mSRef<<" IsTemp="<<a.IsTemp()
+    <<"), mSz="<<mSz<<" mSRef="<<mSRef<<" IsTemp="<<mIsTemp<<endl;
+#endif
+
 if(this == &a) return *this;
 if(a.mSz!=mSz)
@@ -264,8 +311,37 @@
 // NDataBlock a; a = v; ou bien NDataBlock a(10); a = v;
 {
+#ifdef DEBUG_NDATABLOCK
+cout<<"DEBUG_NDataBlock::operator=("<<this<<","<<v<<")"
+    <<" mSz="<<mSz<<" mSRef="<<mSRef<<" IsTemp="<<mIsTemp<<endl;
+#endif
+
 if(mSz==0) throw(SzMismatchError("NDataBlock::operator=v null size\n"));
 T *p=Begin(), *pe=End();
 while (p<pe) *p++ = v;
 return *this;
+}
+
+template <class T>
+T NDataBlock<T>::Sum(size_t i1,size_t n) const
+// Somme des elements de i1 a i1+n-1
+{
+if(i1>=mSz) return 0;
+if(n>mSz) n = mSz; if(n==0) n = mSz-i1;
+T const *p=Begin()+i1, *pe=p+n;
+T val = 0;
+while (p<pe) val += *p++;
+return val;
+}
+
+template <class T>
+T NDataBlock<T>::Product(size_t i1,size_t n) const
+// Produit des elements de i1 a i1+n-1
+{
+if(i1>=mSz) return 0;
+if(n>mSz) n = mSz; if(n==0) n = mSz-i1;
+T const *p=Begin()+i1, *pe=p+n;
+T val = 0;
+while (p<pe) val *= *p++;
+return val;
 }
 
@@ -312,9 +388,10 @@
 
 template <class T>
-NDataBlock<T>& NDataBlock<T>::operator += (NDataBlock<T>& a)
+NDataBlock<T>& NDataBlock<T>::operator += (const NDataBlock<T>& a)
 {
 if(mSz==0 || mSz!=a.mSz)
   throw(SzMismatchError("NDataBlock::operator+=A size mismatch/null"));
-T *p=Begin(), *pe=End(), *pa=a.Begin();
+T *p=Begin(), *pe=End();
+T const * pa=a.Begin();
 while (p<pe) *p++ += *pa++; // ca marche meme si *this=a
 return *this;
@@ -322,9 +399,10 @@
 
 template <class T>
-NDataBlock<T>& NDataBlock<T>::operator -= (NDataBlock<T>& a)
+NDataBlock<T>& NDataBlock<T>::operator -= (const NDataBlock<T>& a)
 {
 if(mSz==0 || mSz!=a.mSz)
   throw(SzMismatchError("NDataBlock::operator-=A size mismatch/null"));
-T *p=Begin(), *pe=End(), *pa=a.Begin();
+T *p=Begin(), *pe=End();
+T const *pa=a.Begin();
 while (p<pe) *p++ -= *pa++; // ca marche meme si *this=a
 return *this;
@@ -332,9 +410,10 @@
 
 template <class T>
-NDataBlock<T>& NDataBlock<T>::operator *= (NDataBlock<T>& a)
+NDataBlock<T>& NDataBlock<T>::operator *= (const NDataBlock<T>& a)
 {
 if(mSz==0 || mSz!=a.mSz)
   throw(SzMismatchError("NDataBlock::operator*=A size mismatch/null"));
-T *p=Begin(), *pe=End(), *pa=a.Begin();
+T *p=Begin(), *pe=End();
+T const *pa=a.Begin();
 while (p<pe) *p++ *= *pa++; // ca marche meme si *this=a
 return *this;
@@ -342,9 +421,10 @@
 
 template <class T>
-NDataBlock<T>& NDataBlock<T>::operator /= (NDataBlock<T>& a)
+NDataBlock<T>& NDataBlock<T>::operator /= (const NDataBlock<T>& a)
 {
 if(mSz==0 || mSz!=a.mSz)
   throw(SzMismatchError("NDataBlock::operator/=A size mismatch/null"));
-T *p=Begin(), *pe=End(), *pa=a.Begin();
+T *p=Begin(), *pe=End();
+T const *pa=a.Begin();
 while (p<pe) *p++ /= *pa++;
 return *this;
@@ -355,28 +435,27 @@
 
 template <class T>
-NDataBlock<T> NDataBlock<T>::Add(T b)
+NDataBlock<T> NDataBlock<T>::Add(T b) const
 // Pour A+b
 {
-NDataBlock<T>& a = *this;
-NDataBlock<T> result(a,false); result.SetTemp(true);
-return result += b;
-}
-
-template <class T>
-NDataBlock<T> NDataBlock<T>::Sub(T b)
+NDataBlock<T> result(*this,false); result.SetTemp(true);
+result += b;
+return result;
+}
+
+template <class T>
+NDataBlock<T> NDataBlock<T>::Sub(T b) const
 // Pour A-b
 {
-NDataBlock<T>& a = *this;
-NDataBlock<T> result(a,false); result.SetTemp(true);
+NDataBlock<T> result(*this,false); result.SetTemp(true);
 return result -= b;
 }
 
 template <class T>
-NDataBlock<T> NDataBlock<T>::SubInv(T b)
+NDataBlock<T> NDataBlock<T>::SubInv(T b) const
 // Pour b-A
 {
-NDataBlock<T>& a = *this;
-NDataBlock<T> result(a,false); result.SetTemp(true);
-T *p=result.Begin(), *pe=result.End(), *pa=a.Begin();
+NDataBlock<T> result(*this,false); result.SetTemp(true);
+T *p=result.Begin(), *pe=result.End();
+T const *pa=this->Begin();
 while(p<pe) {*p++ = b - *pa++;}
 return result;
@@ -384,28 +463,26 @@
 
 template <class T>
-NDataBlock<T> NDataBlock<T>::Mul(T b)
+NDataBlock<T> NDataBlock<T>::Mul(T b) const
 // Pour A*b
 {
-NDataBlock<T>& a = *this;
-NDataBlock<T> result(a,false); result.SetTemp(true);
+NDataBlock<T> result(*this,false); result.SetTemp(true);
 return result *= b;
 }
 
 template <class T>
-NDataBlock<T> NDataBlock<T>::Div(T b)
+NDataBlock<T> NDataBlock<T>::Div(T b) const
 // Pour A/b
 {
-NDataBlock<T>& a = *this;
-NDataBlock<T> result(a,false); result.SetTemp(true);
+NDataBlock<T> result(*this,false); result.SetTemp(true);
 return result /= b;
 }
 
 template <class T>
-NDataBlock<T> NDataBlock<T>::DivInv(T b)
+NDataBlock<T> NDataBlock<T>::DivInv(T b) const
 // Pour b/A
 {
-NDataBlock<T>& a = *this;
-NDataBlock<T> result(a,false); result.SetTemp(true);
-T *p=result.Begin(), *pe=result.End(), *pa=a.Begin();
+NDataBlock<T> result(*this,false); result.SetTemp(true);
+T *p=result.Begin(), *pe=result.End();
+T const *pa = this->Begin();
 while(p<pe) {*p++ = b / *pa++;}
 return result;
@@ -415,69 +492,73 @@
 
 template <class T>
-NDataBlock<T> NDataBlock<T>::Add(NDataBlock<T>& b)
+NDataBlock<T> NDataBlock<T>::Add(const NDataBlock<T>& b) const
 // Pour A+B
 {
-NDataBlock<T>& a = *this;
-if(a.mSz!=b.mSz)
+if(this->mSz!=b.mSz)
   throw(SzMismatchError("NDataBlock operator C=A+B size mismatch/null\n"));
+NDataBlock<T> result; result.SetTemp(true);
 if(b.IsTemp()) {
-  NDataBlock<T> result(b,true); result.SetTemp(true);
-  return result += a;
+  result.Share(b);
+  result += *this;
 } else {
-  NDataBlock<T> result(a,false); result.SetTemp(true);
-  return result += b;
-}
-}
-
-template <class T>
-NDataBlock<T> NDataBlock<T>::Mul(NDataBlock<T>& b)
+  result.Clone(*this);
+  result += b;
+}
+return result;
+}
+
+template <class T>
+NDataBlock<T> NDataBlock<T>::Mul(const NDataBlock<T>& b) const
 // Pour A*B
 {
-NDataBlock<T>& a = *this;
-if(a.mSz!=b.mSz)
+if(this->mSz!=b.mSz)
   throw(SzMismatchError("NDataBlock operator C=A*B size mismatch/null\n"));
+NDataBlock<T> result; result.SetTemp(true);
 if(b.IsTemp()) {
-  NDataBlock<T> result(b,true); result.SetTemp(true);
-  return result *= a;
+  result.Share(b);
+  result *= *this;
 } else {
-  NDataBlock<T> result(a,false); result.SetTemp(true);
-  return result *= b;
-}
-}
-
-template <class T>
-NDataBlock<T> NDataBlock<T>::Sub(NDataBlock<T>& b)
+  result.Clone(*this);
+  result *= b;
+}
+return result;
+}
+
+template <class T>
+NDataBlock<T> NDataBlock<T>::Sub(const NDataBlock<T>& b) const
 // Pour A-B
 {
-NDataBlock<T>& a = *this;
-if(a.mSz!=b.mSz)
+if(this->mSz!=b.mSz)
   throw(SzMismatchError("NDataBlock operator C=A-B size mismatch/null\n"));
+NDataBlock<T> result; result.SetTemp(true);
 if(b.IsTemp()) {
-  NDataBlock<T> result(b,true); result.SetTemp(true);
-  T *p=result.Begin(), *pe=result.End(), *pa=a.Begin();
+  result.Share(b);
+  T *p=result.Begin(), *pe=result.End();
+  T const *pa=this->Begin();
   while(p<pe) {*p = *pa++  - *p; p++;}
-  return result;
 } else {
-  NDataBlock<T> result(a,false); result.SetTemp(true);
-  return result -= b;
-}
-}
-
-template <class T>
-NDataBlock<T> NDataBlock<T>::Div(NDataBlock<T>& b)
+  result.Clone(*this);
+  result -= b;
+}
+return result;
+}
+
+template <class T>
+NDataBlock<T> NDataBlock<T>::Div(const NDataBlock<T>& b) const
 // Pour A/B
 {
-NDataBlock<T>& a = *this;
-if(a.mSz!=b.mSz)
+if(this->mSz!=b.mSz)
   throw(SzMismatchError("NDataBlock operator C=A/B size mismatch/null\n"));
+NDataBlock<T> result; result.SetTemp(true);
 if(b.IsTemp()) {
-  NDataBlock<T> result(b,true); result.SetTemp(true);
-  T *p=result.Begin(), *pe=result.End(), *pa=a.Begin();
+  result.Share(b);
+  T *p=result.Begin(), *pe=result.End();
+  T const *pa=this->Begin();
   while(p<pe) {*p = *pa++  / *p; p++;}
-  return result;
 } else {
-  NDataBlock<T> result(a,false); result.SetTemp(true);
-  return result /= b;
-}
+  result.Clone(*this);
+  result /= b;
+}
+return result;
 }
 
@@ -494,7 +575,6 @@
 #pragma define_template NDataBlock<r_4>
 #pragma define_template NDataBlock<r_8>
-// pas de definitions de cout pour les complex ??
-//#pragma define_template NDataBlock< complex<float> >
-//#pragma define_template NDataBlock< complex<double> >
+#pragma define_template NDataBlock< complex<float> >
+#pragma define_template NDataBlock< complex<double> >
 #endif
 
