diff --git a/histogram/dw_histogram.c b/histogram/dw_histogram.c
index d38c54f0d846ebafde554a27591c875828bb99f6..cf870164810daebda8a18409eabdbe91812ba89b 100644
--- a/histogram/dw_histogram.c
+++ b/histogram/dw_histogram.c
@@ -141,15 +141,20 @@ void AddMatrixObservation(TMatrix X, TMatrixHistogram *h)
  h->sample_size++;
 }
 
-void MatrixPercentile(TMatrix X, PRECISION percentile, TMatrixHistogram *h)
+TMatrix MatrixPercentile(TMatrix X, PRECISION percentile, TMatrixHistogram *h)
 {
  int i, j;
 
- if ((h->rows != RowM(X)) || (h->cols != ColM(X))) dw_Error(SIZE_ERR);
+ if (!X)
+   { if (!(X=CreateMatrix(h->rows,h->cols))) return (TMatrix)NULL; }
+ else
+   if ((h->rows != RowM(X)) || (h->cols != ColM(X))) { dw_Error(SIZE_ERR); return (TMatrix)NULL; }
 
  for (i=h->rows-1; i >= 0; i--)
   for (j=h->cols-1; j >= 0; j--)
    ElementM(X,i,j)=Percentile(percentile,h->low[i][j],h->freq[i][j],ElementM(h->Min,i,j),ElementM(h->Max,i,j),h->intervals,h->sample_size);
+
+ return X;
 }
 
 /*
@@ -278,14 +283,19 @@ void AddVectorObservation(TVector x, TVectorHistogram *h)
  h->sample_size++;
 }
 
-void VectorPercentile(TVector x, PRECISION percentile, TVectorHistogram *h)
+TVector VectorPercentile(TVector x, PRECISION percentile, TVectorHistogram *h)
 {
  int i;
 
- if (h->dim != DimV(x)) dw_Error(SIZE_ERR);
+ if (!x)
+   { if (!(x=CreateVector(h->dim))) return (TVector)NULL; }
+ else
+   if (h->dim != DimV(x)) { dw_Error(SIZE_ERR); return (TVector)NULL; }
 
  for (i=h->dim-1; i >= 0; i--)
   ElementV(x,i)=Percentile(percentile,h->low[i],h->freq[i],ElementV(h->Min,i),ElementV(h->Max,i),h->intervals,h->sample_size);
+
+ return x;
 }
 
 
diff --git a/include/dw_histogram.h b/include/dw_histogram.h
index 2fe0b71e37cad9da4f97ba9a4aeba2427c68f32e..4b75d907bb5c3781de61f658e814b8ae671cc0e0 100644
--- a/include/dw_histogram.h
+++ b/include/dw_histogram.h
@@ -69,7 +69,7 @@ TMatrixHistogram *CreateMatrixHistogram(int rows, int cols, int intervals, int t
 void SetMaxMinMatrixHistogram(TMatrix Min, TMatrix Max, TMatrixHistogram *h);
 void FreeMatrixHistogram(TMatrixHistogram *h);
 void AddMatrixObservation(TMatrix X, TMatrixHistogram *h);
-void MatrixPercentile(TMatrix X, PRECISION percentile, TMatrixHistogram *h);
+TMatrix MatrixPercentile(TMatrix X, PRECISION percentile, TMatrixHistogram *h);
 void MatrixCumulative(TMatrix P, TMatrix Level, TMatrixHistogram *h);
 TMatrix PlotMatrixHistogramAuto(int i, int j, int bins, TMatrixHistogram *h);
 TMatrix PlotMatrixHistogram(int i, int j, PRECISION min, PRECISION max, int bins, TMatrixHistogram *h);
@@ -78,7 +78,7 @@ TVectorHistogram *CreateVectorHistogram(int dim, int intervals, int type);
 void SetMaxMinVectorHistogram(TVector Min, TVector Max, TVectorHistogram *h);
 void FreeVectorHistogram(TVectorHistogram *h);
 void AddVectorObservation(TVector X, TVectorHistogram *h);
-void VectorPercentile(TVector X, PRECISION percentile, TVectorHistogram *h);
+TVector VectorPercentile(TVector X, PRECISION percentile, TVectorHistogram *h);
 void VectorCumulative(TVector p, TVector level, TVectorHistogram *h);
 TMatrix PlotVectorHistogramAuto(int i, int bins, TVectorHistogram *h);
 TMatrix PlotVectorHistogram(int i, PRECISION min, PRECISION max, int bins, TVectorHistogram *h);