diff --git a/CMakeLists.txt b/CMakeLists.txt
index e4b80e6d9a62702af711522d8b45ee4692f0493a..c78032a76b02f94ff002d314f85bcc02dcc710f7 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -32,7 +32,7 @@ execute_process( COMMAND ln -sf ${RooUnfoldHeaders} -t ${CMAKE_CURRENT_BINARY_DI
 
 include(cmake/PlainROOT.cmake)
 
-if(${RooUnfoldTests} AND NOT RooUnfoldTests)  
+if(RooUnfoldTests)  
   include(cmake/testing.cmake)
 endif()
 
diff --git a/src/RooUnfold.h b/src/RooUnfold.h
index b99ada2fac131a70e5a5079db02647befe498ed0..7ddfe15e71bad32f6499e03f60d025689199c75c 100644
--- a/src/RooUnfold.h
+++ b/src/RooUnfold.h
@@ -57,7 +57,7 @@ public:
 
   // Special constructors
 
-  RooUnfoldT<Hist,Hist2D> (const RooUnfoldResponseT<Hist,Hist2D>* res, const Hist* meas, const char* name= 0, const char* title= 0);
+  RooUnfoldT(const RooUnfoldResponseT<Hist,Hist2D>* res, const Hist* meas, const char* name= 0, const char* title= 0);
 
   // Set up an existing object
   
diff --git a/src/RooUnfoldBayes.cxx b/src/RooUnfoldBayes.cxx
index c3fe3e666e9c2aa275de3608fa0427e92b124876..db8dd55df6629e3d5b4fc272dcb571c415838abb 100644
--- a/src/RooUnfoldBayes.cxx
+++ b/src/RooUnfoldBayes.cxx
@@ -51,6 +51,29 @@ RooUnfoldBayesT<Hist,Hist2D>::RooUnfoldBayesT (const RooUnfoldResponseT<Hist,His
   Init();
 }
 
+template<class Hist,class Hist2D>
+void RooUnfoldBayesT<Hist,Hist2D>::SetPriors(const TH1 *priors)
+{
+  if(priors == nullptr){
+    std::cout << "Error SetPriors: no priors provided" << std::endl;
+  }
+  Int_t nbins = priors->GetNbinsX();
+  if( nbins != _nc){
+    std::cout << "Error SetPriors: user's prior size:" << nbins << " different from response matrix size:" << _nc << std::endl;
+  }
+  _P0C.Clear();
+  _P0C.ResizeTo(nbins);
+  for(Int_t i = 1; i< nbins+1; i++){
+    _P0C[i-1] = priors->GetBinContent(i);
+  }
+  _priors = 0;
+  if(this->verbose() >= 2) {
+    std::cout << "Priors set to:" << std::endl;
+    _P0C.Print();
+  }
+}
+
+
 template<class Hist,class Hist2D> RooUnfolding::Algorithm
 RooUnfoldBayesT<Hist,Hist2D>::GetAlgorithm () const
 {
@@ -168,10 +191,19 @@ RooUnfoldBayesT<Hist,Hist2D>::setup() const
   if (this->_dosys)    this->_dnCidPjk.ResizeTo(this->_nc,this->_ne*this->_nc);
 
   // Initial distribution
-  this->_N0C= this->_nCi.Sum();
-  if (this->_N0C!=0.0) {
-    this->_P0C= this->_nCi;
-    this->_P0C *= 1.0/this->_N0C;
+  if(_priors==1){
+    this->_N0C= this->_nCi.Sum();
+    if (this->_N0C!=0.0) {
+      this->_P0C= this->_nCi;
+      this->_P0C *= 1.0/this->_N0C;
+    }
+    if(this->verbose() >=0 ) std::cout << "Using response matrix priors" << std::endl;
+  } else {
+    if(this->verbose() >= 0) std::cout << "Using user's priors" << std::endl;
+  }
+  if(this->verbose()>=1) {
+    std::cout << "Priors:" << std::endl;
+    _P0C.Print();
   }
 }
 
diff --git a/src/RooUnfoldBayes.h b/src/RooUnfoldBayes.h
index 1f0a28e66c2967c13eea974dde85bb184b76e1de..582b9fd597d1b0a69c8b92e96cf1f2b1d7a53078 100644
--- a/src/RooUnfoldBayes.h
+++ b/src/RooUnfoldBayes.h
@@ -40,6 +40,7 @@ public:
   Int_t GetIterations() const;
   Int_t GetSmoothing()  const;
   const TMatrixD& UnfoldingMatrix() const;
+  void SetPriors (const TH1* priors);
 
   virtual RooUnfolding::Algorithm GetAlgorithm() const override;  
   virtual void  SetRegParm (Double_t parm) override;
@@ -71,6 +72,7 @@ protected:
   mutable int _niter;
   mutable int _smoothit;
   mutable int _handleFakes;
+  mutable int _priors =1;
 
   mutable int _nc;              // number of causes  (same as _nt)
   mutable int _ne;              // number of effects (same as _nm)
@@ -92,7 +94,7 @@ protected:
   mutable TMatrixD _dnCidPjk;     // response error propagation matrix (stack j,k into each column)
 
 public:
-  ClassDefOverride (RooUnfoldBayesT, 1) // Bayesian Unfolding
+  ClassDefOverride (RooUnfoldBayesT, 2) // Bayesian Unfolding
 };
 
 
diff --git a/src/RooUnfoldResponse.cxx b/src/RooUnfoldResponse.cxx
index 7539289b1dde4b39e48030d8f1c281c9368bfcf4..1314b9a7954032591d319bc8cf62a22b58f2baaa 100644
--- a/src/RooUnfoldResponse.cxx
+++ b/src/RooUnfoldResponse.cxx
@@ -66,7 +66,7 @@ using std::sqrt;
 template<class Hist, class Hist2D>
 class RooUnfoldFoldingFunction {
 public:
-  RooUnfoldFoldingFunction<Hist,Hist2D> (const RooUnfoldResponseT<Hist,Hist2D>* res, TF1* func, Double_t eps=1e-12, bool verbose=false)
+  RooUnfoldFoldingFunction(const RooUnfoldResponseT<Hist,Hist2D>* res, TF1* func, Double_t eps=1e-12, bool verbose=false)
     : _res(res), _func(func), _eps(eps), _verbose(verbose), _fvals(_res->GetNbinsMeasured()) {
     _ndim= dynamic_cast<TF3*>(_func) ? 3 :
       dynamic_cast<TF2*>(_func) ? 2 : 1;