Logo ROOT  

Reference Guide
 
Loading...
Searching...
No Matches
StandardHistFactoryPlotsWithCategories.C File Reference
Tutorials » RooStats Tutorials

Detailed Description

View in nbviewer Open in SWAN
StandardHistFactoryPlotsWithCategories

This is a standard demo that can be used with any ROOT file prepared in the standard way. You specify:

  • name for input ROOT file
  • name of workspace inside ROOT file that holds model and data
  • name of ModelConfig that specifies details for calculator tools
  • name of dataset

With default parameters the macro will attempt to run the standard hist2workspace example and read the ROOT file that it produces.

The macro will scan through all the categories in a simPdf find the corresponding observable. For each category, it will loop through each of the nuisance parameters and plot

  • the data
  • the nominal model (blue)
  • the +Nsigma (red)
  • the -Nsigma (green)

You can specify how many sigma to vary by changing nSigmaToVary. You can also change the signal rate by changing muVal.

The script produces a lot plots, you can merge them by doing:

gs -q -dNOPAUSE -dBATCH -sDEVICE=pdfwrite -sOutputFile=merged.pdf `ls *pdf`
q
float * q
Definition THbookFile.cxx:89
1) 0x558635269430 RooRealVar:: alpha_syst2 = 0 +/- 1 L(-5 - 5) "alpha_syst2"
2) 0x5586378a28b0 RooRealVar:: alpha_syst3 = 0 +/- 1 L(-5 - 5) "alpha_syst3"
3) 0x558634904e40 RooRealVar:: gamma_stat_channel1_bin_0 = 1 +/- 0.05 L(0 - 1.25) "gamma_stat_channel1_bin_0"
4) 0x558637891360 RooRealVar:: gamma_stat_channel1_bin_1 = 1 +/- 0.1 L(0 - 1.5) "gamma_stat_channel1_bin_1"
check expectedData by category
Is a simultaneous PDF
on type channel1
channelCat==channelCat::channel1
channel1 channel1
[#1] INFO:Plotting -- RooTreeData::plotOn: plotting 234 events out of 234 total events
[#1] INFO:NumericIntegration -- RooRealIntegral::init(channel1_model_Int[obs_x_channel1]) using numeric integrator RooBinIntegrator to calculate Int(obs_x_channel1)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(channel1_model_Int[obs_x_channel1]) using numeric integrator RooBinIntegrator to calculate Int(obs_x_channel1)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(channel1_model_Int[obs_x_channel1]) using numeric integrator RooBinIntegrator to calculate Int(obs_x_channel1)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(channel1_model_Int[obs_x_channel1]) using numeric integrator RooBinIntegrator to calculate Int(obs_x_channel1)
channelCat==channelCat::channel1
channel1 channel1
[#1] INFO:Plotting -- RooTreeData::plotOn: plotting 234 events out of 234 total events
[#1] INFO:NumericIntegration -- RooRealIntegral::init(channel1_model_Int[obs_x_channel1]) using numeric integrator RooBinIntegrator to calculate Int(obs_x_channel1)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(channel1_model_Int[obs_x_channel1]) using numeric integrator RooBinIntegrator to calculate Int(obs_x_channel1)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(channel1_model_Int[obs_x_channel1]) using numeric integrator RooBinIntegrator to calculate Int(obs_x_channel1)
channelCat==channelCat::channel1
channel1 channel1
[#1] INFO:Plotting -- RooTreeData::plotOn: plotting 234 events out of 234 total events
[#1] INFO:NumericIntegration -- RooRealIntegral::init(channel1_model_Int[obs_x_channel1]) using numeric integrator RooBinIntegrator to calculate Int(obs_x_channel1)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(channel1_model_Int[obs_x_channel1]) using numeric integrator RooBinIntegrator to calculate Int(obs_x_channel1)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(channel1_model_Int[obs_x_channel1]) using numeric integrator RooBinIntegrator to calculate Int(obs_x_channel1)
channelCat==channelCat::channel1
channel1 channel1
[#1] INFO:Plotting -- RooTreeData::plotOn: plotting 234 events out of 234 total events
[#1] INFO:NumericIntegration -- RooRealIntegral::init(channel1_model_Int[obs_x_channel1]) using numeric integrator RooBinIntegrator to calculate Int(obs_x_channel1)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(channel1_model_Int[obs_x_channel1]) using numeric integrator RooBinIntegrator to calculate Int(obs_x_channel1)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(channel1_model_Int[obs_x_channel1]) using numeric integrator RooBinIntegrator to calculate Int(obs_x_channel1)
#include "TFile.h"
#include "TROOT.h"
#include "TCanvas.h"
#include "TList.h"
#include "TMath.h"
#include "TSystem.h"
#include "RooWorkspace.h"
#include "RooAbsData.h"
#include "RooRealVar.h"
#include "RooPlot.h"
#include "RooSimultaneous.h"
#include "RooCategory.h"
#include "RooStats/ModelConfig.h"
#include "RooStats/ProfileInspector.h"
using namespace RooFit;
using namespace RooStats;
using std::cout, std::endl;
void StandardHistFactoryPlotsWithCategories(const char *infile = "", const char *workspaceName = "combined",
const char *modelConfigName = "ModelConfig",
const char *dataName = "obsData")
{
double nSigmaToVary = 5.;
double muVal = 0;
bool doFit = false;
// -------------------------------------------------------
// First part is just to access a user-defined file
// or create the standard example file if it doesn't exist
const char *filename = "";
if (!strcmp(infile, "")) {
filename = "results/example_combined_GaussExample_model.root";
bool fileExist = !gSystem->AccessPathName(filename); // note opposite return code
// if file does not exists generate with histfactory
if (!fileExist) {
// Normally this would be run on the command line
cout << "will run standard hist2workspace example" << endl;
gROOT->ProcessLine(".! prepareHistFactory .");
gROOT->ProcessLine(".! hist2workspace config/example.xml");
cout << "\n\n---------------------" << endl;
cout << "Done creating example input" << endl;
cout << "---------------------\n\n" << endl;
}
} else
filename = infile;
// Try to open the file
TFile *file = TFile::Open(filename);
// if input file was specified byt not found, quit
if (!file) {
cout << "StandardRooStatsDemoMacro: Input file " << filename << " is not found" << endl;
return;
}
// -------------------------------------------------------
// Tutorial starts here
// -------------------------------------------------------
// get the workspace out of the file
RooWorkspace *w = (RooWorkspace *)file->Get(workspaceName);
if (!w) {
cout << "workspace not found" << endl;
return;
}
// get the modelConfig out of the file
ModelConfig *mc = (ModelConfig *)w->obj(modelConfigName);
// get the modelConfig out of the file
RooAbsData *data = w->data(dataName);
// make sure ingredients are found
if (!data || !mc) {
w->Print();
cout << "data or ModelConfig was not found" << endl;
return;
}
// -------------------------------------------------------
// now use the profile inspector
RooRealVar *obs = (RooRealVar *)mc->GetObservables()->first();
TList *list = new TList();
RooRealVar *firstPOI = dynamic_cast<RooRealVar *>(mc->GetParametersOfInterest()->first());
firstPOI->setVal(muVal);
// firstPOI->setConstant();
if (doFit) {
mc->GetPdf()->fitTo(*data);
}
// -------------------------------------------------------
mc->GetNuisanceParameters()->Print("v");
int nPlotsMax = 1000;
cout << " check expectedData by category" << endl;
RooDataSet *simData = NULL;
RooSimultaneous *simPdf = NULL;
if (strcmp(mc->GetPdf()->ClassName(), "RooSimultaneous") == 0) {
cout << "Is a simultaneous PDF" << endl;
simPdf = (RooSimultaneous *)(mc->GetPdf());
} else {
cout << "Is not a simultaneous PDF" << endl;
}
if (doFit) {
RooCategory *channelCat = (RooCategory *)(&simPdf->indexCat());
auto const& catName = channelCat->begin()->first;
RooAbsPdf *pdftmp = ((RooSimultaneous *)mc->GetPdf())->getPdf(catName.c_str());
std::unique_ptr<RooArgSet> obstmp{pdftmp->getObservables(*mc->GetObservables())};
obs = ((RooRealVar *)obstmp->first());
RooPlot *frame = obs->frame();
cout << Form("%s==%s::%s", channelCat->GetName(), channelCat->GetName(), catName.c_str()) << endl;
cout << catName << " " << channelCat->getLabel() << endl;
data->plotOn(frame, MarkerSize(1),
Cut(Form("%s==%s::%s", channelCat->GetName(), channelCat->GetName(), catName.c_str())),
DataError(RooAbsData::None));
Double_t normCount =
data->sumEntries(Form("%s==%s::%s", channelCat->GetName(), channelCat->GetName(), catName.c_str()));
pdftmp->plotOn(frame, LineWidth(2.), Normalization(normCount, RooAbsReal::NumEvent));
frame->Draw();
cout << "expected events = " << mc->GetPdf()->expectedEvents(*data->get()) << endl;
return;
}
int nPlots = 0;
if (!simPdf) {
for (auto *var : static_range_cast<RooRealVar *>(*mc->GetNuisanceParameters())) {
RooPlot *frame = obs->frame();
frame->SetYTitle(var->GetName());
data->plotOn(frame, MarkerSize(1));
var->setVal(0);
mc->GetPdf()->plotOn(frame, LineWidth(1.));
var->setVal(1);
mc->GetPdf()->plotOn(frame, LineColor(kRed), LineStyle(kDashed), LineWidth(1));
var->setVal(-1);
mc->GetPdf()->plotOn(frame, LineColor(kGreen), LineStyle(kDashed), LineWidth(1));
list->Add(frame);
var->setVal(0);
}
} else {
RooCategory *channelCat = (RooCategory *)(&simPdf->indexCat());
for (auto const& tt : *channelCat) {
if (nPlots == nPlotsMax) {
break;
}
auto const& catName = tt.first;
cout << "on type " << catName << " " << endl;
// Get pdf associated with state from simpdf
RooAbsPdf *pdftmp = simPdf->getPdf(catName.c_str());
// Generate observables defined by the pdf associated with this state
std::unique_ptr<RooArgSet> obstmp{pdftmp->getObservables(*mc->GetObservables())};
// obstmp->Print();
obs = ((RooRealVar *)obstmp->first());
for (auto *var : static_range_cast<RooRealVar*>(*mc->GetNuisanceParameters())) {
if (nPlots == nPlotsMax) break;
TCanvas *c2 = new TCanvas("c2");
RooPlot *frame = obs->frame();
frame->SetName(Form("frame%d", nPlots));
frame->SetYTitle(var->GetName());
cout << Form("%s==%s::%s", channelCat->GetName(), channelCat->GetName(), catName.c_str()) << endl;
cout << catName << " " << channelCat->getLabel() << endl;
data->plotOn(frame, MarkerSize(1),
Cut(Form("%s==%s::%s", channelCat->GetName(), channelCat->GetName(), catName.c_str())),
DataError(RooAbsData::None));
Double_t normCount =
data->sumEntries(Form("%s==%s::%s", channelCat->GetName(), channelCat->GetName(), catName.c_str()));
if (strcmp(var->GetName(), "Lumi") == 0) {
cout << "working on lumi" << endl;
var->setVal(w->var("nominalLumi")->getVal());
var->Print();
} else {
var->setVal(0);
}
// w->allVars().Print("v");
// mc->GetNuisanceParameters()->Print("v");
// pdftmp->plotOn(frame,LineWidth(2.));
// mc->GetPdf()->plotOn(frame,LineWidth(2.),Slice(*channelCat,catName.c_str()),ProjWData(*data));
// pdftmp->plotOn(frame,LineWidth(2.),Slice(*channelCat,catName.c_str()),ProjWData(*data));
normCount = pdftmp->expectedEvents(*obs);
pdftmp->plotOn(frame, LineWidth(2.), Normalization(normCount, RooAbsReal::NumEvent));
if (strcmp(var->GetName(), "Lumi") == 0) {
cout << "working on lumi" << endl;
var->setVal(w->var("nominalLumi")->getVal() + 0.05);
var->Print();
} else {
var->setVal(nSigmaToVary);
}
// pdftmp->plotOn(frame,LineColor(kRed),LineStyle(kDashed),LineWidth(2));
// mc->GetPdf()->plotOn(frame,LineColor(kRed),LineStyle(kDashed),LineWidth(2.),Slice(*channelCat,catName.c_str()),ProjWData(*data));
// pdftmp->plotOn(frame,LineColor(kRed),LineStyle(kDashed),LineWidth(2.),Slice(*channelCat,catName.c_str()),ProjWData(*data));
normCount = pdftmp->expectedEvents(*obs);
pdftmp->plotOn(frame, LineWidth(2.), LineColor(kRed), LineStyle(kDashed),
Normalization(normCount, RooAbsReal::NumEvent));
if (strcmp(var->GetName(), "Lumi") == 0) {
cout << "working on lumi" << endl;
var->setVal(w->var("nominalLumi")->getVal() - 0.05);
var->Print();
} else {
var->setVal(-nSigmaToVary);
}
// pdftmp->plotOn(frame,LineColor(kGreen),LineStyle(kDashed),LineWidth(2));
// mc->GetPdf()->plotOn(frame,LineColor(kGreen),LineStyle(kDashed),LineWidth(2),Slice(*channelCat,catName.c_str()),ProjWData(*data));
// pdftmp->plotOn(frame,LineColor(kGreen),LineStyle(kDashed),LineWidth(2),Slice(*channelCat,catName.c_str()),ProjWData(*data));
normCount = pdftmp->expectedEvents(*obs);
pdftmp->plotOn(frame, LineWidth(2.), LineColor(kGreen), LineStyle(kDashed),
Normalization(normCount, RooAbsReal::NumEvent));
// set them back to normal
if (strcmp(var->GetName(), "Lumi") == 0) {
cout << "working on lumi" << endl;
var->setVal(w->var("nominalLumi")->getVal());
var->Print();
} else {
var->setVal(0);
}
list->Add(frame);
// quit making plots
++nPlots;
frame->Draw();
c2->SaveAs(Form("%s_%s_%s.pdf", catName.c_str(), obs->GetName(), var->GetName()));
delete c2;
}
}
}
// -------------------------------------------------------
// now make plots
TCanvas *c1 = new TCanvas("c1", "ProfileInspectorDemo", 800, 200);
if (list->GetSize() > 4) {
double n = list->GetSize();
int nx = (int)sqrt(n);
int ny = TMath::CeilNint(n / nx);
nx = TMath::CeilNint(sqrt(n));
c1->Divide(ny, nx);
} else
c1->Divide(list->GetSize());
for (int i = 0; i < list->GetSize(); ++i) {
c1->cd(i + 1);
list->At(i)->Draw();
}
}
static_range_cast
ROOT::RRangeCast< T, false, Range_t > static_range_cast(Range_t &&coll)
Definition RooAbsCollection.h:45
Double_t
double Double_t
Definition RtypesCore.h:59
kRed
@ kRed
Definition Rtypes.h:66
kGreen
@ kGreen
Definition Rtypes.h:66
kDashed
@ kDashed
Definition TAttLine.h:48
data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Definition TGWin32VirtualXProxy.cxx:104
filename
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char filename
Definition TGWin32VirtualXProxy.cxx:232
gROOT
#define gROOT
Definition TROOT.h:406
Form
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition TString.cxx:2489
gSystem
R__EXTERN TSystem * gSystem
Definition TSystem.h:561
RooAbsArg::getObservables
RooFit::OwningPtr< RooArgSet > getObservables(const RooArgSet &set, bool valueOnly=true) const
Given a set of possible observables, return the observables that this PDF depends on.
Definition RooAbsArg.cxx:698
RooAbsCategory::getLabel
const char * getLabel() const
Retrieve current label. Use getCurrentLabel() for more clarity.
Definition RooAbsCategory.h:135
RooAbsCategory::begin
std::map< std::string, value_type >::const_iterator begin() const
Iterator for category state names. Points to pairs of index and name.
Definition RooAbsCategory.h:91
RooAbsCollection::first
RooAbsArg * first() const
Definition RooAbsCollection.h:291
RooAbsCollection::Print
void Print(Option_t *options=nullptr) const override
This method must be overridden when a class wants to print itself.
Definition RooAbsCollection.h:315
RooAbsData
Abstract base class for binned and unbinned datasets.
Definition RooAbsData.h:57
RooAbsPdf
Abstract interface for all probability density functions.
Definition RooAbsPdf.h:40
RooAbsPdf::expectedEvents
virtual double expectedEvents(const RooArgSet *nset) const
Return expected number of events to be used in calculation of extended likelihood.
Definition RooAbsPdf.cxx:2398
RooAbsPdf::plotOn
RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1={}, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={}, const RooCmdArg &arg6={}, const RooCmdArg &arg7={}, const RooCmdArg &arg8={}, const RooCmdArg &arg9={}, const RooCmdArg &arg10={}) const override
Helper calling plotOn(RooPlot*, RooLinkedList&) const.
Definition RooAbsPdf.h:124
RooAbsPdf::fitTo
RooFit::OwningPtr< RooFitResult > fitTo(RooAbsData &data, CmdArgs_t const &... cmdArgs)
Fit PDF to given dataset.
Definition RooAbsPdf.h:157
RooAbsRealLValue::frame
RooPlot * frame(const RooCmdArg &arg1, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={}, const RooCmdArg &arg6={}, const RooCmdArg &arg7={}, const RooCmdArg &arg8={}) const
Create a new RooPlot on the heap with a drawing frame initialized for this object,...
Definition RooAbsRealLValue.cxx:216
RooCategory
Object to represent discrete states.
Definition RooCategory.h:28
RooDataSet
Container class to hold unbinned data.
Definition RooDataSet.h:33
RooPlot
Plot frame and a container for graphics objects within that frame.
Definition RooPlot.h:45
RooPlot::Print
void Print(Option_t *options=nullptr) const override
This method must be overridden when a class wants to print itself.
Definition RooPlot.h:138
RooPlot::Draw
void Draw(Option_t *options=nullptr) override
Draw this plot and all of the elements it contains.
Definition RooPlot.cxx:637
RooPlot::SetYTitle
void SetYTitle(const char *title)
Definition RooPlot.cxx:1336
RooPlot::SetName
Use the constructor that doesn t take the name and and then call SetName() and SetTitle() on the RooPlot.")
RooRealVar
Variable that can be changed from the outside.
Definition RooRealVar.h:37
RooRealVar::setVal
void setVal(double value) override
Set value of variable to 'value'.
Definition RooRealVar.cxx:242
RooSimultaneous
Facilitates simultaneous fitting of multiple PDFs to subsets of a given dataset.
Definition RooSimultaneous.h:41
RooSimultaneous::getPdf
RooAbsPdf * getPdf(RooStringView catName) const
Return the p.d.f associated with the given index category name.
Definition RooSimultaneous.cxx:368
RooSimultaneous::indexCat
const RooAbsCategoryLValue & indexCat() const
Definition RooSimultaneous.h:91
RooStats::ModelConfig
ModelConfig is a simple class that holds configuration information specifying how a model should be u...
Definition ModelConfig.h:35
RooStats::ModelConfig::GetParametersOfInterest
const RooArgSet * GetParametersOfInterest() const
get RooArgSet containing the parameter of interest (return nullptr if not existing)
Definition ModelConfig.h:274
RooStats::ModelConfig::GetNuisanceParameters
const RooArgSet * GetNuisanceParameters() const
get RooArgSet containing the nuisance parameters (return nullptr if not existing)
Definition ModelConfig.h:277
RooStats::ModelConfig::GetObservables
const RooArgSet * GetObservables() const
get RooArgSet for observables (return nullptr if not existing)
Definition ModelConfig.h:286
RooStats::ModelConfig::GetPdf
RooAbsPdf * GetPdf() const
get model PDF (return nullptr if pdf has not been specified or does not exist)
Definition ModelConfig.h:271
RooWorkspace
Persistable container for RooFit projects.
Definition RooWorkspace.h:43
TCanvas
The Canvas class.
Definition TCanvas.h:23
TCollection::GetSize
virtual Int_t GetSize() const
Return the capacity of the collection, i.e.
Definition TCollection.h:184
TDirectoryFile::Get
TObject * Get(const char *namecycle) override
Return pointer to object identified by namecycle.
Definition TDirectoryFile.cxx:937
TFile
A ROOT file is an on-disk file, usually with extension .root, that stores objects in a file-system-li...
Definition TFile.h:53
TFile::Open
static TFile * Open(const char *name, Option_t *option="", const char *ftitle="", Int_t compress=ROOT::RCompressionSetting::EDefaults::kUseCompiledDefault, Int_t netopt=0)
Create / open a file.
Definition TFile.cxx:4089
TList
A doubly linked list.
Definition TList.h:38
TList::Add
void Add(TObject *obj) override
Definition TList.h:81
TList::At
TObject * At(Int_t idx) const override
Returns the object at position idx. Returns 0 if idx is out of range.
Definition TList.cxx:355
TNamed::GetName
const char * GetName() const override
Returns name of object.
Definition TNamed.h:47
TObject::ClassName
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition TObject.cxx:213
TObject::Draw
virtual void Draw(Option_t *option="")
Default Draw method for all objects.
Definition TObject.cxx:280
TSystem::AccessPathName
virtual Bool_t AccessPathName(const char *path, EAccessMode mode=kFileExists)
Returns FALSE if one can access a file using the specified access mode.
Definition TSystem.cxx:1296
RooFit::DataError
RooCmdArg DataError(Int_t)
Definition RooGlobalFunc.cxx:470
RooFit::LineWidth
RooCmdArg LineWidth(Width_t width)
Definition RooGlobalFunc.cxx:248
RooFit::Normalization
RooCmdArg Normalization(double scaleFactor)
Definition RooGlobalFunc.cxx:195
RooFit::MarkerSize
RooCmdArg MarkerSize(Size_t size)
Definition RooGlobalFunc.cxx:348
RooFit::Cut
RooCmdArg Cut(const char *cutSpec)
Definition RooGlobalFunc.cxx:320
RooFit::LineColor
RooCmdArg LineColor(Color_t color)
Definition RooGlobalFunc.cxx:215
RooFit::LineStyle
RooCmdArg LineStyle(Style_t style)
Definition RooGlobalFunc.cxx:240
c1
return c1
Definition legend1.C:41
n
const Int_t n
Definition legend1.C:16
c2
return c2
Definition legend2.C:14
ROOT::Math::sqrt
VecExpr< UnaryOp< Sqrt< T >, VecExpr< A, T, D >, T >, T, D > sqrt(const VecExpr< A, T, D > &rhs)
Definition UnaryOperators.h:281
RooFit
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Definition JSONIO.h:26
RooStats
Namespace for the RooStats classes.
Definition Asimov.h:19
TMath::CeilNint
Int_t CeilNint(Double_t x)
Returns the nearest integer of TMath::Ceil(x).
Definition TMath.h:674
tt
auto * tt
Definition textangle.C:16
Author
Kyle Cranmer

Definition in file StandardHistFactoryPlotsWithCategories.C.