diff --git a/Code/BasicFilters/otbStreamingStatisticsVectorImageFilter.txx b/Code/BasicFilters/otbStreamingStatisticsVectorImageFilter.txx
index 2a21c80f1eddd5832431f386b4456f3724803ea7..a04b6f37dd80adb774c57feb8d1e275465ec5b13 100644
--- a/Code/BasicFilters/otbStreamingStatisticsVectorImageFilter.txx
+++ b/Code/BasicFilters/otbStreamingStatisticsVectorImageFilter.txx
@@ -290,7 +290,6 @@ namespace otb {
     // sum of squares 
     for( i = 0; i < numberOfThreads; i++)
       {
-	std::cout<<i<<std::endl;
 	count += m_Count[i];
 	/** TODO
 	 *  To modify using + method operator. If we use it now -> exceptionmacro (no GetClassName...)
@@ -351,7 +350,7 @@ namespace otb {
     tempTranspose = pixelSumMatrix.GetTranspose();
     covMatrixTemp = pixelSumMatrix*tempTranspose;
     /** TODO
-     * A modifier en utilisant l'opérateur + de la méthode. Pour le moment probleme avec exceptionmacro (pas de GetClassName...)
+     *  To modify using - method operator. If we use it now -> exceptionmacro (no GetClassName...)
      *covMatrix -= covMatrixTemp;
      **/
     if( (covMatrix.Rows() != covMatrixTemp.Rows()) || (covMatrix.Cols() != covMatrixTemp.Cols()))
@@ -430,9 +429,8 @@ namespace otb {
 	++it;
 	progress.CompletedPixel();
 	pixelTransposeVector = pixelVector.GetTranspose();
-      
-	/** TODO
-	 * A modifier en utilisant l'opérateur + de la méthode. Pour le moment probleme avec exceptionmacro (pas de GetClassName...)
+        /** TODO
+         *  To modify using + method operator. If we use it now -> exceptionmacro (no GetClassName...)
 	 * m_XX[threadId]+=pixelVector*pixelTransposeVector;
 	 **/
 	tempMatrix = pixelVector*pixelTransposeVector;
diff --git a/Code/Fusion/otbBayesianFusionFilter.h b/Code/Fusion/otbBayesianFusionFilter.h
index 0773a1777a34377ccffdd57bcf326a3b15e524ad..1652c5496b1ef04a8bd4ccefb5d0918648727d0b 100755
--- a/Code/Fusion/otbBayesianFusionFilter.h
+++ b/Code/Fusion/otbBayesianFusionFilter.h
@@ -245,15 +245,18 @@ public:
   itkSetMacro(S, float);
   /** Give the S coefficient. */
   itkGetConstReferenceMacro(S, float);
-
-
+  
 protected:
   BayesianFusionFilter();
   virtual ~BayesianFusionFilter();
-  /** Initialize some accumulators before the threads run. */
+  /** Check if internal statistics need to be computed, and do so */
   void BeforeThreadedGenerateData ();
- 
-                        
+  /** Compute internal statistics required for fusion */
+  void ComputeInternalStatistics(void);
+  /** Call the superclass implementation and set the StatisticsHaveBeenGenerated
+   * flag to false */
+  virtual void Modified(void);
+   
 private:
   /** Ponderation declaration*/
   float m_Lambda;
@@ -266,7 +269,8 @@ private:
   MatrixType m_Beta;
   /** Optimisation matrix */
   MatrixType m_Vcondopt;
-  
+  /** True if internal statistics have been generated */
+  bool m_StatisticsHaveBeenGenerated;
 };
 
 } // end namespace otb
diff --git a/Code/Fusion/otbBayesianFusionFilter.txx b/Code/Fusion/otbBayesianFusionFilter.txx
index 335328bceb926221b74fd6fefed18fd4a773156e..beeecd68f48f2fc0c396195a4b485d3fb560b849 100755
--- a/Code/Fusion/otbBayesianFusionFilter.txx
+++ b/Code/Fusion/otbBayesianFusionFilter.txx
@@ -40,6 +40,7 @@ namespace otb
   {
     m_Lambda = 0.9999;
     m_S = 1;
+     m_StatisticsHaveBeenGenerated = false;
   }
   
   
@@ -55,7 +56,20 @@ namespace otb
   {
     
   }
-  
+  template <class TInputMultiSpectralImage, 
+  class TInputMultiSpectralInterpImage, 
+  class TInputPanchroImage,  
+  class TOutputImage>
+      void
+      BayesianFusionFilter<TInputMultiSpectralImage, 
+      TInputMultiSpectralInterpImage,  
+      TInputPanchroImage, 
+      TOutputImage>
+  ::Modified()
+      {
+       Superclass::Modified();
+        m_StatisticsHaveBeenGenerated = false; 
+      }
   
   template <class TInputMultiSpectralImage, 
 	    class TInputMultiSpectralInterpImage, 
@@ -68,7 +82,25 @@ namespace otb
 		       TOutputImage>
   ::BeforeThreadedGenerateData ()
   {	
-    OutputImageRegionType msiRequestedRegion = this->GetMultiSpectInterp()->GetRequestedRegion();
+    if(!m_StatisticsHaveBeenGenerated)
+    {
+     this->ComputeInternalStatistics();
+     m_StatisticsHaveBeenGenerated = true;
+    }
+  }
+  
+  template <class TInputMultiSpectralImage, 
+  class TInputMultiSpectralInterpImage, 
+  class TInputPanchroImage,  
+  class TOutputImage>
+      void
+      BayesianFusionFilter<TInputMultiSpectralImage, 
+      TInputMultiSpectralInterpImage,  
+      TInputPanchroImage, 
+      TOutputImage>
+  ::ComputeInternalStatistics()
+        {	
+      OutputImageRegionType msiRequestedRegion = this->GetMultiSpectInterp()->GetRequestedRegion();
     OutputImageRegionType msRequestedRegion = this->GetMultiSpect()->GetRequestedRegion();
     OutputImageRegionType panchroRequestedRegion = this->GetPanchro()->GetRequestedRegion();
 
@@ -95,8 +127,10 @@ namespace otb
 
     covComputefilter->SetInput(multiSpecInterp);
     covComputefilter->Update();
-
+    
+    
     MatrixType m_CovarianceMatrix = covComputefilter->GetCovariance();
+    otbMsgDebugMacro(<<"Covariance: "<<m_CovarianceMatrix);
 
     m_CovarianceInvMatrix = m_CovarianceMatrix.GetInverse();
     /** Beta computation : Regression model coefficient */
@@ -120,8 +154,10 @@ namespace otb
     msTransposePan->SetSecondInput( caster->GetOutput() );
     
     msTransposeMs->Update();
+    otbMsgDebugMacro(<<"MsTMs: "<<msTransposeMs->GetResultOutput()->Get());
     msTransposePan->Update();
-
+    otbMsgDebugMacro(<<"MsTPan: "<<msTransposePan->GetResultOutput()->Get()); 
+    
     MatrixType temp;
     temp = msTransposeMs->GetResultOutput()->Get().GetInverse();
     m_Beta = temp*msTransposePan->GetResultOutput()->Get();
@@ -132,6 +168,7 @@ namespace otb
     panTransposePan->SetFirstInput(caster->GetOutput());
     panTransposePan->SetSecondInput(caster->GetOutput());
     panTransposePan->Update();
+    otbMsgDebugMacro(<<"PanTPan: "<<msTransposePan->GetResultOutput()->Get()); 
     MatrixType S, tempS, tempS2;
     S = panTransposePan->GetResultOutput()->Get();
     tempS = msTransposePan->GetResultOutput()->Get().GetTranspose();
@@ -282,6 +319,8 @@ namespace otb
    panchro->SetRequestedRegion(panchroRequestedRegion);
    panchro->PropagateRequestedRegion();
    panchro->UpdateOutputData(); 
+   
+    std::cout<<"End of BeforeThreaderGenerateData."<<std::endl;
   } 
 } // end namespace otb