ROOT logo
/*****************************************************************************

 *****************************************************************************/

//////////////////////////////////////////////////////////////////////////////
// 
// BEGIN_HTML
// RooJeffreysPrior 
// END_HTML
//


#include "RooFit.h"

#include "Riostream.h"
#include "Riostream.h"
#include <math.h>

#include "RooJeffreysPrior.h"
#include "RooAbsReal.h"
#include "RooAbsPdf.h"
#include "RooErrorHandler.h"
#include "RooArgSet.h"
#include "RooMsgService.h"
#include "RooFitResult.h"
#include "TMatrixDSym.h"
#include "RooDataHist.h"
#include "RooFitResult.h"
#include "RooNumIntConfig.h"
#include "RooRealVar.h"

ClassImp(RooJeffreysPrior)
;

using namespace RooFit;

//_____________________________________________________________________________
RooJeffreysPrior::RooJeffreysPrior(const char* name, const char* title, 
			     RooAbsPdf& nominal,
			     const RooArgList& paramSet,
			     const RooArgList& obsSet) :
  RooAbsPdf(name, title),
  _nominal("nominal","nominal",this,nominal,kFALSE,kFALSE),
  _obsSet("!obsSet","obs-side variation",this,kFALSE,kFALSE),
  _paramSet("!paramSet","high-side variation",this)
{
  //_obsSet("!obsSet","obs-side variation",this),
  _obsIter = _obsSet.createIterator() ;
  _paramIter = _paramSet.createIterator() ;


  TIterator* inputIter1 = obsSet.createIterator() ;
  RooAbsArg* comp ;
  while((comp = (RooAbsArg*)inputIter1->Next())) {
    if (!dynamic_cast<RooAbsReal*>(comp)) {
      coutE(InputArguments) << "RooJeffreysPrior::ctor(" << GetName() << ") ERROR: component " << comp->GetName() 
			    << " in first list is not of type RooAbsReal" << endl ;
      RooErrorHandler::softAbort() ;
    }
    _obsSet.add(*comp) ;
    //    if (takeOwnership) {
    //      _ownedList.addOwned(*comp) ;
    //    }
  }
  delete inputIter1 ;



  TIterator* inputIter3 = paramSet.createIterator() ;
  while((comp = (RooAbsArg*)inputIter3->Next())) {
    if (!dynamic_cast<RooAbsReal*>(comp)) {
      coutE(InputArguments) << "RooJeffreysPrior::ctor(" << GetName() << ") ERROR: component " << comp->GetName() 
			    << " in first list is not of type RooAbsReal" << endl ;
      RooErrorHandler::softAbort() ;
    }
    _paramSet.add(*comp) ;
    //    if (takeOwnership) {
    //      _ownedList.addOwned(*comp) ;
    //    }
  }
  delete inputIter3 ;


  // use a different integrator by default.
  if(paramSet.getSize()==1)
    this->specialIntegratorConfig(kTRUE)->method1D().setLabel("RooAdaptiveGaussKronrodIntegrator1D")  ;

}



//_____________________________________________________________________________
RooJeffreysPrior::RooJeffreysPrior(const RooJeffreysPrior& other, const char* name) :
  RooAbsPdf(other, name), 
  _nominal("!nominal",this,other._nominal),
  _obsSet("!obsSet",this,other._obsSet),
  _paramSet("!paramSet",this,other._paramSet)
{
  // Copy constructor
  _obsIter = _obsSet.createIterator() ;
  _paramIter = _paramSet.createIterator() ;

  // Member _ownedList is intentionally not copy-constructed -- ownership is not transferred
}

//_____________________________________________________________________________
RooJeffreysPrior::RooJeffreysPrior() 
{
  // Default constructor
  _obsIter = NULL;
  _paramIter = NULL;

}



//_____________________________________________________________________________
RooJeffreysPrior::~RooJeffreysPrior() 
{
  // Destructor

  if (_obsIter) delete _obsIter ;
  if (_paramIter) delete _paramIter ;
}




//_____________________________________________________________________________
Double_t RooJeffreysPrior::evaluate() const 
{
  // Calculate and return current value of self
  RooFit::MsgLevel msglevel = RooMsgService::instance().globalKillBelow();
  RooMsgService::instance().setGlobalKillBelow(RooFit::WARNING);
  // create Asimov dataset
  //  _paramSet.Print("v");
  //  return sqrt(_paramSet.getRealValue("mu"));
  *(_nominal.arg().getVariables()) = _paramSet;
  /*
  cout << "_______________"<<endl;
  _paramSet.Print("v");
  _nominal->getVariables()->Print("v");
  cout << "_______________"<<endl;
  */
  RooDataHist* data = ((RooAbsPdf&)(_nominal.arg())).generateBinned(_obsSet,ExpectedData());
  //  data->Print("v");
  //RooFitResult* res = _nominal->fitTo(*data, Save(),PrintLevel(-1),Minos(kFALSE),SumW2Error(kTRUE));
  RooFitResult* res = ((RooAbsPdf&)(_nominal.arg())).fitTo(*data, Save(),PrintLevel(-1),Minos(kFALSE),SumW2Error(kFALSE));
  TMatrixDSym cov = res->covarianceMatrix();
  cov.Invert();
  double ret =  sqrt(cov.Determinant());
  
  /*
    // for 1 parameter can avoid making TMatrix etc.
    // but number of params may be > 1 with others held constant
  if(_paramSet.getSize()==1){
    RooRealVar* var = (RooRealVar*) _paramSet.first();
    // also, the _paramSet proxy one does not pick up a different value
    cout << "eval at "<< ret << " " << 1/(var->getError()) << endl; 
    // need to get the actual variable instance out of the pdf like below
    var = (RooRealVar*) _nominal->getVariables()->find(var->GetName());
    cout << "eval at "<< ret << " " << 1/(var->getError()) << endl; 
  }
  */

  //  res->Print();
  delete data;
  delete res;
  RooMsgService::instance().setGlobalKillBelow(msglevel);

  //  cout << "eval at "<< ret << endl; 
  //  _paramSet.Print("v");
  return ret;

}

//_____________________________________________________________________________
Int_t RooJeffreysPrior::getAnalyticalIntegral(RooArgSet& /*allVars*/, RooArgSet& /*analVars*/, const char* /*rangeName*/) const 
{
  //  if (matchArgs(allVars,analVars,x)) return 1 ;
  //  if (matchArgs(allVars,analVars,mean)) return 2 ;
  //  return 1;
  return 0 ;
}



//_____________________________________________________________________________
Double_t RooJeffreysPrior::analyticalIntegral(Int_t code, const char* /*rangeName*/) const 
{
  assert(code==1 );
  //cout << "evaluating analytic integral" << endl;
  return 1.;
}


 RooJeffreysPrior.cxx:1
 RooJeffreysPrior.cxx:2
 RooJeffreysPrior.cxx:3
 RooJeffreysPrior.cxx:4
 RooJeffreysPrior.cxx:5
 RooJeffreysPrior.cxx:6
 RooJeffreysPrior.cxx:7
 RooJeffreysPrior.cxx:8
 RooJeffreysPrior.cxx:9
 RooJeffreysPrior.cxx:10
 RooJeffreysPrior.cxx:11
 RooJeffreysPrior.cxx:12
 RooJeffreysPrior.cxx:13
 RooJeffreysPrior.cxx:14
 RooJeffreysPrior.cxx:15
 RooJeffreysPrior.cxx:16
 RooJeffreysPrior.cxx:17
 RooJeffreysPrior.cxx:18
 RooJeffreysPrior.cxx:19
 RooJeffreysPrior.cxx:20
 RooJeffreysPrior.cxx:21
 RooJeffreysPrior.cxx:22
 RooJeffreysPrior.cxx:23
 RooJeffreysPrior.cxx:24
 RooJeffreysPrior.cxx:25
 RooJeffreysPrior.cxx:26
 RooJeffreysPrior.cxx:27
 RooJeffreysPrior.cxx:28
 RooJeffreysPrior.cxx:29
 RooJeffreysPrior.cxx:30
 RooJeffreysPrior.cxx:31
 RooJeffreysPrior.cxx:32
 RooJeffreysPrior.cxx:33
 RooJeffreysPrior.cxx:34
 RooJeffreysPrior.cxx:35
 RooJeffreysPrior.cxx:36
 RooJeffreysPrior.cxx:37
 RooJeffreysPrior.cxx:38
 RooJeffreysPrior.cxx:39
 RooJeffreysPrior.cxx:40
 RooJeffreysPrior.cxx:41
 RooJeffreysPrior.cxx:42
 RooJeffreysPrior.cxx:43
 RooJeffreysPrior.cxx:44
 RooJeffreysPrior.cxx:45
 RooJeffreysPrior.cxx:46
 RooJeffreysPrior.cxx:47
 RooJeffreysPrior.cxx:48
 RooJeffreysPrior.cxx:49
 RooJeffreysPrior.cxx:50
 RooJeffreysPrior.cxx:51
 RooJeffreysPrior.cxx:52
 RooJeffreysPrior.cxx:53
 RooJeffreysPrior.cxx:54
 RooJeffreysPrior.cxx:55
 RooJeffreysPrior.cxx:56
 RooJeffreysPrior.cxx:57
 RooJeffreysPrior.cxx:58
 RooJeffreysPrior.cxx:59
 RooJeffreysPrior.cxx:60
 RooJeffreysPrior.cxx:61
 RooJeffreysPrior.cxx:62
 RooJeffreysPrior.cxx:63
 RooJeffreysPrior.cxx:64
 RooJeffreysPrior.cxx:65
 RooJeffreysPrior.cxx:66
 RooJeffreysPrior.cxx:67
 RooJeffreysPrior.cxx:68
 RooJeffreysPrior.cxx:69
 RooJeffreysPrior.cxx:70
 RooJeffreysPrior.cxx:71
 RooJeffreysPrior.cxx:72
 RooJeffreysPrior.cxx:73
 RooJeffreysPrior.cxx:74
 RooJeffreysPrior.cxx:75
 RooJeffreysPrior.cxx:76
 RooJeffreysPrior.cxx:77
 RooJeffreysPrior.cxx:78
 RooJeffreysPrior.cxx:79
 RooJeffreysPrior.cxx:80
 RooJeffreysPrior.cxx:81
 RooJeffreysPrior.cxx:82
 RooJeffreysPrior.cxx:83
 RooJeffreysPrior.cxx:84
 RooJeffreysPrior.cxx:85
 RooJeffreysPrior.cxx:86
 RooJeffreysPrior.cxx:87
 RooJeffreysPrior.cxx:88
 RooJeffreysPrior.cxx:89
 RooJeffreysPrior.cxx:90
 RooJeffreysPrior.cxx:91
 RooJeffreysPrior.cxx:92
 RooJeffreysPrior.cxx:93
 RooJeffreysPrior.cxx:94
 RooJeffreysPrior.cxx:95
 RooJeffreysPrior.cxx:96
 RooJeffreysPrior.cxx:97
 RooJeffreysPrior.cxx:98
 RooJeffreysPrior.cxx:99
 RooJeffreysPrior.cxx:100
 RooJeffreysPrior.cxx:101
 RooJeffreysPrior.cxx:102
 RooJeffreysPrior.cxx:103
 RooJeffreysPrior.cxx:104
 RooJeffreysPrior.cxx:105
 RooJeffreysPrior.cxx:106
 RooJeffreysPrior.cxx:107
 RooJeffreysPrior.cxx:108
 RooJeffreysPrior.cxx:109
 RooJeffreysPrior.cxx:110
 RooJeffreysPrior.cxx:111
 RooJeffreysPrior.cxx:112
 RooJeffreysPrior.cxx:113
 RooJeffreysPrior.cxx:114
 RooJeffreysPrior.cxx:115
 RooJeffreysPrior.cxx:116
 RooJeffreysPrior.cxx:117
 RooJeffreysPrior.cxx:118
 RooJeffreysPrior.cxx:119
 RooJeffreysPrior.cxx:120
 RooJeffreysPrior.cxx:121
 RooJeffreysPrior.cxx:122
 RooJeffreysPrior.cxx:123
 RooJeffreysPrior.cxx:124
 RooJeffreysPrior.cxx:125
 RooJeffreysPrior.cxx:126
 RooJeffreysPrior.cxx:127
 RooJeffreysPrior.cxx:128
 RooJeffreysPrior.cxx:129
 RooJeffreysPrior.cxx:130
 RooJeffreysPrior.cxx:131
 RooJeffreysPrior.cxx:132
 RooJeffreysPrior.cxx:133
 RooJeffreysPrior.cxx:134
 RooJeffreysPrior.cxx:135
 RooJeffreysPrior.cxx:136
 RooJeffreysPrior.cxx:137
 RooJeffreysPrior.cxx:138
 RooJeffreysPrior.cxx:139
 RooJeffreysPrior.cxx:140
 RooJeffreysPrior.cxx:141
 RooJeffreysPrior.cxx:142
 RooJeffreysPrior.cxx:143
 RooJeffreysPrior.cxx:144
 RooJeffreysPrior.cxx:145
 RooJeffreysPrior.cxx:146
 RooJeffreysPrior.cxx:147
 RooJeffreysPrior.cxx:148
 RooJeffreysPrior.cxx:149
 RooJeffreysPrior.cxx:150
 RooJeffreysPrior.cxx:151
 RooJeffreysPrior.cxx:152
 RooJeffreysPrior.cxx:153
 RooJeffreysPrior.cxx:154
 RooJeffreysPrior.cxx:155
 RooJeffreysPrior.cxx:156
 RooJeffreysPrior.cxx:157
 RooJeffreysPrior.cxx:158
 RooJeffreysPrior.cxx:159
 RooJeffreysPrior.cxx:160
 RooJeffreysPrior.cxx:161
 RooJeffreysPrior.cxx:162
 RooJeffreysPrior.cxx:163
 RooJeffreysPrior.cxx:164
 RooJeffreysPrior.cxx:165
 RooJeffreysPrior.cxx:166
 RooJeffreysPrior.cxx:167
 RooJeffreysPrior.cxx:168
 RooJeffreysPrior.cxx:169
 RooJeffreysPrior.cxx:170
 RooJeffreysPrior.cxx:171
 RooJeffreysPrior.cxx:172
 RooJeffreysPrior.cxx:173
 RooJeffreysPrior.cxx:174
 RooJeffreysPrior.cxx:175
 RooJeffreysPrior.cxx:176
 RooJeffreysPrior.cxx:177
 RooJeffreysPrior.cxx:178
 RooJeffreysPrior.cxx:179
 RooJeffreysPrior.cxx:180
 RooJeffreysPrior.cxx:181
 RooJeffreysPrior.cxx:182
 RooJeffreysPrior.cxx:183
 RooJeffreysPrior.cxx:184
 RooJeffreysPrior.cxx:185
 RooJeffreysPrior.cxx:186
 RooJeffreysPrior.cxx:187
 RooJeffreysPrior.cxx:188
 RooJeffreysPrior.cxx:189
 RooJeffreysPrior.cxx:190
 RooJeffreysPrior.cxx:191
 RooJeffreysPrior.cxx:192
 RooJeffreysPrior.cxx:193
 RooJeffreysPrior.cxx:194
 RooJeffreysPrior.cxx:195
 RooJeffreysPrior.cxx:196