From $ROOTSYS/tutorials/math/testUnfold1.C

// Author: Stefan Schmitt
// DESY, 14.10.2008

//  Version 17.0, updated for using the classes TUnfoldDensity, TUnfoldBinning
//
//  History:
//    Version 16.1, parallel to changes in TUnfold
//    Version 16.0, parallel to changes in TUnfold
//    Version 15, with automated L-curve scan
//    Version 14, with changes in TUnfoldSys.cxx
//    Version 13, include test of systematic errors
//    Version 12, catch error when defining the input
//    Version 11,  print chi**2 and number of degrees of freedom
//    Version 10,  with bug-fix in TUnfold.cxx
//    Version 9,  with bug-fix in TUnfold.cxx and TUnfold.h
//    Version 8,  with bug-fix in TUnfold.cxx and TUnfold.h
//    Version 7,  with bug-fix in TUnfold.cxx and TUnfold.h
//    Version 6a, fix problem with dynamic array allocation under windows
//    Version 6, bug-fixes in TUnfold.C
//    Version 5, replace main() by testUnfold1()
//    Version 4, with bug-fix in TUnfold.C
//    Version 3, with bug-fix in TUnfold.C
//    Version 2, with changed ScanLcurve() arguments
//    Version 1, remove L curve analysis, use ScanLcurve() method instead
//    Version 0, L curve analysis included here

/*
  This file is part of TUnfold.

  TUnfold is free software: you can redistribute it and/or modify
  it under the terms of the GNU General Public License as published by
  the Free Software Foundation, either version 3 of the License, or
  (at your option) any later version.

  TUnfold is distributed in the hope that it will be useful,
  but WITHOUT ANY WARRANTY; without even the implied warranty of
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  GNU General Public License for more details.

  You should have received a copy of the GNU General Public License
  along with TUnfold.  If not, see <http://www.gnu.org/licenses/>.
*/

#include <TError.h>
#include <TMath.h>
#include <TCanvas.h>
#include <TRandom3.h>
#include <TFitter.h>
#include <TF1.h>
#include <TStyle.h>
#include <TVector.h>
#include <TGraph.h>

#include "TUnfoldDensity.h"

// #define VERBOSE_LCURVE_SCAN

using namespace std;

///////////////////////////////////////////////////////////////////////
//
//  Test program for the classes TUnfold and related
//
//  (1) Generate Monte Carlo and Data events
//      The events consist of
//        signal
//        background
//
//      The signal is a resonance. It is generated with a Breit-Wigner,
//      smeared by a Gaussian
//
//  (2) Unfold the data. The result is:
//      The background level
//      The shape of the resonance, corrected for detector effects
//
//      Systematic errors from the MC shape variation are included
//      and propagated to the result
//
//  (3) fit the unfolded distribution, including the correlation matrix
//
//  (4) save six plots to a file testUnfold1.ps
//        1  2  3
//        4  5  6
//      1: 2d-plot of the matrix decsribing the migrations
//      2: generator-level distributions
//             blue: unfolded data, total errors
//             green: unfolded data, statistical errors
//             red: generated data
//             black: fit to green data points
//      3: detector level distributions
//             blue: unfoldede data, folded back through the matrix
//             black: Monte Carlo (with wrong peal position)
//             blue: data
//      4: global correlation coefficients
//      5: chi**2 as a function of log(tau)
//           the star indicates the final choice of tau
//      6: the L curve
//
///////////////////////////////////////////////////////////////////////

TRandom *rnd=0;

TH2 *gHistInvEMatrix;

TVirtualFitter *gFitter=0;

void chisquare_corr(Int_t &npar, Double_t * /*gin */, Double_t &f, Double_t *u, Int_t /*flag */) {
  //  Minimization function for H1s using a Chisquare method
  //  only one-dimensional histograms are supported
  //  Corelated errors are taken from an external inverse covariance matrix
  //  stored in a 2-dimensional histogram
  Double_t x;
  TH1 *hfit = (TH1*)gFitter->GetObjectFit();
  TF1 *f1   = (TF1*)gFitter->GetUserFunc();

  f1->InitArgs(&x,u);
  npar = f1->GetNpar();
  f = 0;

  Int_t npfit = 0;
  Int_t nPoints=hfit->GetNbinsX();
  Double_t *df=new Double_t[nPoints];
  for (Int_t i=0;i<nPoints;i++) {
    x     = hfit->GetBinCenter(i+1);
    TF1::RejectPoint(kFALSE);
    df[i] = f1->EvalPar(&x,u)-hfit->GetBinContent(i+1);
    if (TF1::RejectedPoint()) df[i]=0.0;
    else npfit++;
  }
  for (Int_t i=0;i<nPoints;i++) {
    for (Int_t j=0;j<nPoints;j++) {
      f += df[i]*df[j]*gHistInvEMatrix->GetBinContent(i+1,j+1);
    }
  }
  delete[] df;
  f1->SetNumberFitPoints(npfit);
}

Double_t bw_func(Double_t *x,Double_t *par) {
  Double_t dm=x[0]-par[1];
  return par[0]/(dm*dm+par[2]*par[2]);
}


// generate an event
// output:
//  negative mass: background event
//  positive mass: signal event
Double_t GenerateEvent(Double_t bgr, // relative fraction of background
                       Double_t mass, // peak position
                       Double_t gamma) // peak width
{
  Double_t t;
  if(rnd->Rndm()>bgr) {
    // generate signal event
    // with positive mass
    do {
      do {
        t=rnd->Rndm();
      } while(t>=1.0);
      t=TMath::Tan((t-0.5)*TMath::Pi())*gamma+mass;
    } while(t<=0.0);
    return t;
  } else {
    // generate background event
    // generate events following a power-law distribution
    //   f(E) = K * TMath::power((E0+E),N0)
    static Double_t const E0=2.4;
    static Double_t const N0=2.9;
    do {
      do {
        t=rnd->Rndm();
      } while(t>=1.0);
      // the mass is returned negative
      // In our example a convenient way to indicate it is a background event.
      t= -(TMath::Power(1.-t,1./(1.-N0))-1.0)*E0;
    } while(t>=0.0);
    return t;
  }
}

// smear the event to detector level
// input:
//   mass on generator level (mTrue>0 !)
// output:
//   mass on detector level
Double_t DetectorEvent(Double_t mTrue) {
  // smear by double-gaussian
  static Double_t frac=0.1;
  static Double_t wideBias=0.03;
  static Double_t wideSigma=0.5;
  static Double_t smallBias=0.0;
  static Double_t smallSigma=0.1;
  if(rnd->Rndm()>frac) {
    return rnd->Gaus(mTrue+smallBias,smallSigma);
  } else {
    return rnd->Gaus(mTrue+wideBias,wideSigma);
  }
}

int testUnfold1()
{
  // switch on histogram errors
  TH1::SetDefaultSumw2();

  // show fit result
  gStyle->SetOptFit(1111);

  // random generator
  rnd=new TRandom3();

  // data and MC luminosity, cross-section
  Double_t const luminosityData=100000;
  Double_t const luminosityMC=1000000;
  Double_t const crossSection=1.0;

  Int_t const nDet=250;
  Int_t const nGen=100;
  Double_t const xminDet=0.0;
  Double_t const xmaxDet=10.0;
  Double_t const xminGen=0.0;
  Double_t const xmaxGen=10.0;

  //============================================
  // generate MC distribution
  //
  TH1D *histMgenMC=new TH1D("MgenMC",";mass(gen)",nGen,xminGen,xmaxGen);
  TH1D *histMdetMC=new TH1D("MdetMC",";mass(det)",nDet,xminDet,xmaxDet);
  TH2D *histMdetGenMC=new TH2D("MdetgenMC",";mass(det);mass(gen)",
                               nDet,xminDet,xmaxDet,nGen,xminGen,xmaxGen);
  Int_t neventMC=rnd->Poisson(luminosityMC*crossSection);
  for(Int_t i=0;i<neventMC;i++) {
    Double_t mGen=GenerateEvent(0.3, // relative fraction of background
                                4.0, // peak position in MC
                                0.2); // peak width in MC
    Double_t mDet=DetectorEvent(TMath::Abs(mGen));
    // the generated mass is negative for background
    // and positive for signal
    // so it will be filled in the underflow bin
    // this is very convenient for the unfolding:
    // the unfolded result will contain the number of background
    // events in the underflow bin

    // generated MC distribution (for comparison only)
    histMgenMC->Fill(mGen,luminosityData/luminosityMC);
    // reconstructed MC distribution (for comparison only)
    histMdetMC->Fill(mDet,luminosityData/luminosityMC);

    // matrix describing how the generator input migrates to the
    // reconstructed level. Unfolding input.
    // NOTE on underflow/overflow bins:
    //  (1) the detector level under/overflow bins are used for
    //       normalisation ("efficiency" correction)
    //       in our toy example, these bins are populated from tails
    //       of the initial MC distribution.
    //  (2) the generator level underflow/overflow bins are
    //       unfolded. In this example:
    //       underflow bin: background events reconstructed in the detector
    //       overflow bin: signal events generated at masses > xmaxDet
    // for the unfolded result these bins will be filled
    //  -> the background normalisation will be contained in the underflow bin
    histMdetGenMC->Fill(mDet,mGen,luminosityData/luminosityMC);
  }

  //============================================
  // generate alternative MC
  // this will be used to derive a systematic error due to MC
  // parameter uncertainties
  TH2D *histMdetGenSysMC=new TH2D("MdetgenSysMC",";mass(det);mass(gen)",
                                  nDet,xminDet,xmaxDet,nGen,xminGen,xmaxGen);
  neventMC=rnd->Poisson(luminosityMC*crossSection);
  for(Int_t i=0;i<neventMC;i++) {
    Double_t mGen=GenerateEvent
       (0.5, // relative fraction of background
        3.6, // peak position in MC with systematic shift
        0.15); // peak width in MC
    Double_t mDet=DetectorEvent(TMath::Abs(mGen));
    histMdetGenSysMC->Fill(mDet,mGen,luminosityData/luminosityMC);
  }

  //============================================
  // generate data distribution
  //
  TH1D *histMgenData=new TH1D("MgenData",";mass(gen)",nGen,xminGen,xmaxGen);
  TH1D *histMdetData=new TH1D("MdetData",";mass(det)",nDet,xminDet,xmaxDet);
  Int_t neventData=rnd->Poisson(luminosityData*crossSection);
  for(Int_t i=0;i<neventData;i++) {
    Double_t mGen=GenerateEvent(0.4, // relative fraction of background
                                3.8, // peak position in data
                                0.15); // peak width in data
    Double_t mDet=DetectorEvent(TMath::Abs(mGen));
    // generated data mass for comparison plots
    // for real data, we do not have this histogram
    histMgenData->Fill(mGen);

    // reconstructed mass, unfolding input
    histMdetData->Fill(mDet);
  }

  //=========================================================================
  // divide by bin withd to get density distributions
  TH1D *histDensityGenData=new TH1D("DensityGenData",";mass(gen)",
                                    nGen,xminGen,xmaxGen);
  TH1D *histDensityGenMC=new TH1D("DensityGenMC",";mass(gen)",
                                    nGen,xminGen,xmaxGen);
  for(Int_t i=1;i<=nGen;i++) {
     histDensityGenData->SetBinContent(i,histMgenData->GetBinContent(i)/
                                       histMgenData->GetBinWidth(i));
     histDensityGenMC->SetBinContent(i,histMgenMC->GetBinContent(i)/
                                       histMgenMC->GetBinWidth(i));
  }

  //=========================================================================
  // set up the unfolding
  // define migration matrix
  TUnfoldDensity unfold(histMdetGenMC,TUnfold::kHistMapOutputVert);

  // define input and bias scame
  // do not use the bias, because MC peak may be at the wrong place
  // watch out for error codes returned by the SetInput method
  // errors larger or equal 10000 are fatal:
  // the data points specified as input are not sufficient to constrain the
  // unfolding process
  if(unfold.SetInput(histMdetData)>=10000) {
    std::cout<<"Unfolding result may be wrong\n";
  }

  //========================================================================
  // the unfolding is done here
  //
  // scan L curve and find best point
  Int_t nScan=30;
  // use automatic L-curve scan: start with taumin=taumax=0.0
  Double_t tauMin=0.0;
  Double_t tauMax=0.0;
  Int_t iBest;
  TSpline *logTauX,*logTauY;
  TGraph *lCurve;

  // if required, report Info messages (for debugging the L-curve scan)
#ifdef VERBOSE_LCURVE_SCAN
  Int_t oldinfo=gErrorIgnoreLevel;
  gErrorIgnoreLevel=kInfo;
#endif
  // this method scans the parameter tau and finds the kink in the L curve
  // finally, the unfolding is done for the best choice of tau
  iBest=unfold.ScanLcurve(nScan,tauMin,tauMax,&lCurve,&logTauX,&logTauY);

  // if required, switch to previous log-level
#ifdef VERBOSE_LCURVE_SCAN
  gErrorIgnoreLevel=oldinfo;
#endif

  //==========================================================================
  // define a correlated systematic error
  // for example, assume there is a 10% correlated error for all reconstructed
  // masses larger than 7
  Double_t SYS_ERROR1_MSTART=6;
  Double_t SYS_ERROR1_SIZE=0.1;
  TH2D *histMdetGenSys1=new TH2D("Mdetgensys1",";mass(det);mass(gen)",
                                 nDet,xminDet,xmaxDet,nGen,xminGen,xmaxGen);
  for(Int_t i=0;i<=nDet+1;i++) {
     if(histMdetData->GetBinCenter(i)>=SYS_ERROR1_MSTART) {
        for(Int_t j=0;j<=nGen+1;j++) {
           histMdetGenSys1->SetBinContent(i,j,SYS_ERROR1_SIZE);
        }
     }
  }
  unfold.AddSysError(histMdetGenSysMC,"SYSERROR_MC",TUnfold::kHistMapOutputVert,
                     TUnfoldSys::kSysErrModeMatrix);
  unfold.AddSysError(histMdetGenSys1,"SYSERROR1",TUnfold::kHistMapOutputVert,
                     TUnfoldSys::kSysErrModeRelative);

  //==========================================================================
  // print some results
  //
  std::cout<<"tau="<<unfold.GetTau()<<"\n";
  std::cout<<"chi**2="<<unfold.GetChi2A()<<"+"<<unfold.GetChi2L()
           <<" / "<<unfold.GetNdf()<<"\n";
  std::cout<<"chi**2(sys)="<<unfold.GetChi2Sys()<<"\n";


  //==========================================================================
  // create graphs with one point to visualize the best choice of tau
  //
  Double_t t[1],x[1],y[1];
  logTauX->GetKnot(iBest,t[0],x[0]);
  logTauY->GetKnot(iBest,t[0],y[0]);
  TGraph *bestLcurve=new TGraph(1,x,y);
  TGraph *bestLogTauLogChi2=new TGraph(1,t,x);

  //==========================================================================
  // retreive results into histograms

  // get unfolded distribution
  TH1 *histMunfold=unfold.GetOutput("Unfolded");

  // get unfolding result, folded back
  TH1 *histMdetFold=unfold.GetFoldedOutput("FoldedBack");

  // get error matrix (input distribution [stat] errors only)
  // TH2D *histEmatData=unfold.GetEmatrix("EmatData");

  // get total error matrix:
  //   migration matrix uncorrelated and correlated systematic errors
  //   added in quadrature to the data statistical errors
  TH2 *histEmatTotal=unfold.GetEmatrixTotal("EmatTotal");

  // create data histogram with the total errors
  TH1D *histTotalError=
     new TH1D("TotalError",";mass(gen)",nGen,xminGen,xmaxGen);
  for(Int_t bin=1;bin<=nGen;bin++) {
    histTotalError->SetBinContent(bin,histMunfold->GetBinContent(bin));
    histTotalError->SetBinError
       (bin,TMath::Sqrt(histEmatTotal->GetBinContent(bin,bin)));
  }

  // get global correlation coefficients
  // for this calculation one has to specify whether the
  // underflow/overflow bins are included or not
  // default: include all bins
  // here: exclude underflow and overflow bins
  TH1 *histRhoi=unfold.GetRhoItotal("rho_I",
                                    0, // use default title
                                    0, // all distributions
                                    "*[UO]", // discard underflow and overflow bins on all axes
                                    kTRUE, // use original binning
                                    &gHistInvEMatrix // store inverse of error matrix
                                    );

  //======================================================================
  // fit Breit-Wigner shape to unfolded data, using the full error matrix
  // here we use a "user" chi**2 function to take into account
  // the full covariance matrix

  gFitter=TVirtualFitter::Fitter(histMunfold);
  gFitter->SetFCN(chisquare_corr);

  TF1 *bw=new TF1("bw",bw_func,xminGen,xmaxGen,3);
  bw->SetParameter(0,1000.);
  bw->SetParameter(1,3.8);
  bw->SetParameter(2,0.2);

  // for (wrong!) fitting without correlations, drop the option "U"
  // here.
  histMunfold->Fit(bw,"UE");

  //=====================================================================
  // plot some histograms
  TCanvas output;
  output.Divide(3,2);

  // Show the matrix which connects input and output
  // There are overflow bins at the bottom, not shown in the plot
  // These contain the background shape.
  // The overflow bins to the left and right contain
  // events which are not reconstructed. These are necessary for proper MC
  // normalisation
  output.cd(1);
  histMdetGenMC->Draw("BOX");

  // draw generator-level distribution:
  //   data (red) [for real data this is not available]
  //   MC input (black) [with completely wrong peak position and shape]
  //   unfolded data (blue)
  output.cd(2);
  histTotalError->SetLineColor(kBlue);
  histTotalError->Draw("E");
  histMunfold->SetLineColor(kGreen);
  histMunfold->Draw("SAME E1");
  histDensityGenData->SetLineColor(kRed);
  histDensityGenData->Draw("SAME");
  histDensityGenMC->Draw("SAME HIST");

  // show detector level distributions
  //    data (red)
  //    MC (black) [with completely wrong peak position and shape]
  //    unfolded data (blue)
  output.cd(3);
  histMdetFold->SetLineColor(kBlue);
  histMdetFold->Draw();
  histMdetMC->Draw("SAME HIST");

  TH1 *histInput=unfold.GetInput("Minput",";mass(det)");

  histInput->SetLineColor(kRed);
  histInput->Draw("SAME");

  // show correlation coefficients
  output.cd(4);
  histRhoi->Draw();

  // show tau as a function of chi**2
  output.cd(5);
  logTauX->Draw();
  bestLogTauLogChi2->SetMarkerColor(kRed);
  bestLogTauLogChi2->Draw("*");

  // show the L curve
  output.cd(6);
  lCurve->Draw("AL");
  bestLcurve->SetMarkerColor(kRed);
  bestLcurve->Draw("*");

  output.SaveAs("testUnfold1.ps");

  return 0;
}
 testUnfold1.C:1
 testUnfold1.C:2
 testUnfold1.C:3
 testUnfold1.C:4
 testUnfold1.C:5
 testUnfold1.C:6
 testUnfold1.C:7
 testUnfold1.C:8
 testUnfold1.C:9
 testUnfold1.C:10
 testUnfold1.C:11
 testUnfold1.C:12
 testUnfold1.C:13
 testUnfold1.C:14
 testUnfold1.C:15
 testUnfold1.C:16
 testUnfold1.C:17
 testUnfold1.C:18
 testUnfold1.C:19
 testUnfold1.C:20
 testUnfold1.C:21
 testUnfold1.C:22
 testUnfold1.C:23
 testUnfold1.C:24
 testUnfold1.C:25
 testUnfold1.C:26
 testUnfold1.C:27
 testUnfold1.C:28
 testUnfold1.C:29
 testUnfold1.C:30
 testUnfold1.C:31
 testUnfold1.C:32
 testUnfold1.C:33
 testUnfold1.C:34
 testUnfold1.C:35
 testUnfold1.C:36
 testUnfold1.C:37
 testUnfold1.C:38
 testUnfold1.C:39
 testUnfold1.C:40
 testUnfold1.C:41
 testUnfold1.C:42
 testUnfold1.C:43
 testUnfold1.C:44
 testUnfold1.C:45
 testUnfold1.C:46
 testUnfold1.C:47
 testUnfold1.C:48
 testUnfold1.C:49
 testUnfold1.C:50
 testUnfold1.C:51
 testUnfold1.C:52
 testUnfold1.C:53
 testUnfold1.C:54
 testUnfold1.C:55
 testUnfold1.C:56
 testUnfold1.C:57
 testUnfold1.C:58
 testUnfold1.C:59
 testUnfold1.C:60
 testUnfold1.C:61
 testUnfold1.C:62
 testUnfold1.C:63
 testUnfold1.C:64
 testUnfold1.C:65
 testUnfold1.C:66
 testUnfold1.C:67
 testUnfold1.C:68
 testUnfold1.C:69
 testUnfold1.C:70
 testUnfold1.C:71
 testUnfold1.C:72
 testUnfold1.C:73
 testUnfold1.C:74
 testUnfold1.C:75
 testUnfold1.C:76
 testUnfold1.C:77
 testUnfold1.C:78
 testUnfold1.C:79
 testUnfold1.C:80
 testUnfold1.C:81
 testUnfold1.C:82
 testUnfold1.C:83
 testUnfold1.C:84
 testUnfold1.C:85
 testUnfold1.C:86
 testUnfold1.C:87
 testUnfold1.C:88
 testUnfold1.C:89
 testUnfold1.C:90
 testUnfold1.C:91
 testUnfold1.C:92
 testUnfold1.C:93
 testUnfold1.C:94
 testUnfold1.C:95
 testUnfold1.C:96
 testUnfold1.C:97
 testUnfold1.C:98
 testUnfold1.C:99
 testUnfold1.C:100
 testUnfold1.C:101
 testUnfold1.C:102
 testUnfold1.C:103
 testUnfold1.C:104
 testUnfold1.C:105
 testUnfold1.C:106
 testUnfold1.C:107
 testUnfold1.C:108
 testUnfold1.C:109
 testUnfold1.C:110
 testUnfold1.C:111
 testUnfold1.C:112
 testUnfold1.C:113
 testUnfold1.C:114
 testUnfold1.C:115
 testUnfold1.C:116
 testUnfold1.C:117
 testUnfold1.C:118
 testUnfold1.C:119
 testUnfold1.C:120
 testUnfold1.C:121
 testUnfold1.C:122
 testUnfold1.C:123
 testUnfold1.C:124
 testUnfold1.C:125
 testUnfold1.C:126
 testUnfold1.C:127
 testUnfold1.C:128
 testUnfold1.C:129
 testUnfold1.C:130
 testUnfold1.C:131
 testUnfold1.C:132
 testUnfold1.C:133
 testUnfold1.C:134
 testUnfold1.C:135
 testUnfold1.C:136
 testUnfold1.C:137
 testUnfold1.C:138
 testUnfold1.C:139
 testUnfold1.C:140
 testUnfold1.C:141
 testUnfold1.C:142
 testUnfold1.C:143
 testUnfold1.C:144
 testUnfold1.C:145
 testUnfold1.C:146
 testUnfold1.C:147
 testUnfold1.C:148
 testUnfold1.C:149
 testUnfold1.C:150
 testUnfold1.C:151
 testUnfold1.C:152
 testUnfold1.C:153
 testUnfold1.C:154
 testUnfold1.C:155
 testUnfold1.C:156
 testUnfold1.C:157
 testUnfold1.C:158
 testUnfold1.C:159
 testUnfold1.C:160
 testUnfold1.C:161
 testUnfold1.C:162
 testUnfold1.C:163
 testUnfold1.C:164
 testUnfold1.C:165
 testUnfold1.C:166
 testUnfold1.C:167
 testUnfold1.C:168
 testUnfold1.C:169
 testUnfold1.C:170
 testUnfold1.C:171
 testUnfold1.C:172
 testUnfold1.C:173
 testUnfold1.C:174
 testUnfold1.C:175
 testUnfold1.C:176
 testUnfold1.C:177
 testUnfold1.C:178
 testUnfold1.C:179
 testUnfold1.C:180
 testUnfold1.C:181
 testUnfold1.C:182
 testUnfold1.C:183
 testUnfold1.C:184
 testUnfold1.C:185
 testUnfold1.C:186
 testUnfold1.C:187
 testUnfold1.C:188
 testUnfold1.C:189
 testUnfold1.C:190
 testUnfold1.C:191
 testUnfold1.C:192
 testUnfold1.C:193
 testUnfold1.C:194
 testUnfold1.C:195
 testUnfold1.C:196
 testUnfold1.C:197
 testUnfold1.C:198
 testUnfold1.C:199
 testUnfold1.C:200
 testUnfold1.C:201
 testUnfold1.C:202
 testUnfold1.C:203
 testUnfold1.C:204
 testUnfold1.C:205
 testUnfold1.C:206
 testUnfold1.C:207
 testUnfold1.C:208
 testUnfold1.C:209
 testUnfold1.C:210
 testUnfold1.C:211
 testUnfold1.C:212
 testUnfold1.C:213
 testUnfold1.C:214
 testUnfold1.C:215
 testUnfold1.C:216
 testUnfold1.C:217
 testUnfold1.C:218
 testUnfold1.C:219
 testUnfold1.C:220
 testUnfold1.C:221
 testUnfold1.C:222
 testUnfold1.C:223
 testUnfold1.C:224
 testUnfold1.C:225
 testUnfold1.C:226
 testUnfold1.C:227
 testUnfold1.C:228
 testUnfold1.C:229
 testUnfold1.C:230
 testUnfold1.C:231
 testUnfold1.C:232
 testUnfold1.C:233
 testUnfold1.C:234
 testUnfold1.C:235
 testUnfold1.C:236
 testUnfold1.C:237
 testUnfold1.C:238
 testUnfold1.C:239
 testUnfold1.C:240
 testUnfold1.C:241
 testUnfold1.C:242
 testUnfold1.C:243
 testUnfold1.C:244
 testUnfold1.C:245
 testUnfold1.C:246
 testUnfold1.C:247
 testUnfold1.C:248
 testUnfold1.C:249
 testUnfold1.C:250
 testUnfold1.C:251
 testUnfold1.C:252
 testUnfold1.C:253
 testUnfold1.C:254
 testUnfold1.C:255
 testUnfold1.C:256
 testUnfold1.C:257
 testUnfold1.C:258
 testUnfold1.C:259
 testUnfold1.C:260
 testUnfold1.C:261
 testUnfold1.C:262
 testUnfold1.C:263
 testUnfold1.C:264
 testUnfold1.C:265
 testUnfold1.C:266
 testUnfold1.C:267
 testUnfold1.C:268
 testUnfold1.C:269
 testUnfold1.C:270
 testUnfold1.C:271
 testUnfold1.C:272
 testUnfold1.C:273
 testUnfold1.C:274
 testUnfold1.C:275
 testUnfold1.C:276
 testUnfold1.C:277
 testUnfold1.C:278
 testUnfold1.C:279
 testUnfold1.C:280
 testUnfold1.C:281
 testUnfold1.C:282
 testUnfold1.C:283
 testUnfold1.C:284
 testUnfold1.C:285
 testUnfold1.C:286
 testUnfold1.C:287
 testUnfold1.C:288
 testUnfold1.C:289
 testUnfold1.C:290
 testUnfold1.C:291
 testUnfold1.C:292
 testUnfold1.C:293
 testUnfold1.C:294
 testUnfold1.C:295
 testUnfold1.C:296
 testUnfold1.C:297
 testUnfold1.C:298
 testUnfold1.C:299
 testUnfold1.C:300
 testUnfold1.C:301
 testUnfold1.C:302
 testUnfold1.C:303
 testUnfold1.C:304
 testUnfold1.C:305
 testUnfold1.C:306
 testUnfold1.C:307
 testUnfold1.C:308
 testUnfold1.C:309
 testUnfold1.C:310
 testUnfold1.C:311
 testUnfold1.C:312
 testUnfold1.C:313
 testUnfold1.C:314
 testUnfold1.C:315
 testUnfold1.C:316
 testUnfold1.C:317
 testUnfold1.C:318
 testUnfold1.C:319
 testUnfold1.C:320
 testUnfold1.C:321
 testUnfold1.C:322
 testUnfold1.C:323
 testUnfold1.C:324
 testUnfold1.C:325
 testUnfold1.C:326
 testUnfold1.C:327
 testUnfold1.C:328
 testUnfold1.C:329
 testUnfold1.C:330
 testUnfold1.C:331
 testUnfold1.C:332
 testUnfold1.C:333
 testUnfold1.C:334
 testUnfold1.C:335
 testUnfold1.C:336
 testUnfold1.C:337
 testUnfold1.C:338
 testUnfold1.C:339
 testUnfold1.C:340
 testUnfold1.C:341
 testUnfold1.C:342
 testUnfold1.C:343
 testUnfold1.C:344
 testUnfold1.C:345
 testUnfold1.C:346
 testUnfold1.C:347
 testUnfold1.C:348
 testUnfold1.C:349
 testUnfold1.C:350
 testUnfold1.C:351
 testUnfold1.C:352
 testUnfold1.C:353
 testUnfold1.C:354
 testUnfold1.C:355
 testUnfold1.C:356
 testUnfold1.C:357
 testUnfold1.C:358
 testUnfold1.C:359
 testUnfold1.C:360
 testUnfold1.C:361
 testUnfold1.C:362
 testUnfold1.C:363
 testUnfold1.C:364
 testUnfold1.C:365
 testUnfold1.C:366
 testUnfold1.C:367
 testUnfold1.C:368
 testUnfold1.C:369
 testUnfold1.C:370
 testUnfold1.C:371
 testUnfold1.C:372
 testUnfold1.C:373
 testUnfold1.C:374
 testUnfold1.C:375
 testUnfold1.C:376
 testUnfold1.C:377
 testUnfold1.C:378
 testUnfold1.C:379
 testUnfold1.C:380
 testUnfold1.C:381
 testUnfold1.C:382
 testUnfold1.C:383
 testUnfold1.C:384
 testUnfold1.C:385
 testUnfold1.C:386
 testUnfold1.C:387
 testUnfold1.C:388
 testUnfold1.C:389
 testUnfold1.C:390
 testUnfold1.C:391
 testUnfold1.C:392
 testUnfold1.C:393
 testUnfold1.C:394
 testUnfold1.C:395
 testUnfold1.C:396
 testUnfold1.C:397
 testUnfold1.C:398
 testUnfold1.C:399
 testUnfold1.C:400
 testUnfold1.C:401
 testUnfold1.C:402
 testUnfold1.C:403
 testUnfold1.C:404
 testUnfold1.C:405
 testUnfold1.C:406
 testUnfold1.C:407
 testUnfold1.C:408
 testUnfold1.C:409
 testUnfold1.C:410
 testUnfold1.C:411
 testUnfold1.C:412
 testUnfold1.C:413
 testUnfold1.C:414
 testUnfold1.C:415
 testUnfold1.C:416
 testUnfold1.C:417
 testUnfold1.C:418
 testUnfold1.C:419
 testUnfold1.C:420
 testUnfold1.C:421
 testUnfold1.C:422
 testUnfold1.C:423
 testUnfold1.C:424
 testUnfold1.C:425
 testUnfold1.C:426
 testUnfold1.C:427
 testUnfold1.C:428
 testUnfold1.C:429
 testUnfold1.C:430
 testUnfold1.C:431
 testUnfold1.C:432
 testUnfold1.C:433
 testUnfold1.C:434
 testUnfold1.C:435
 testUnfold1.C:436
 testUnfold1.C:437
 testUnfold1.C:438
 testUnfold1.C:439
 testUnfold1.C:440
 testUnfold1.C:441
 testUnfold1.C:442
 testUnfold1.C:443
 testUnfold1.C:444
 testUnfold1.C:445
 testUnfold1.C:446
 testUnfold1.C:447
 testUnfold1.C:448
 testUnfold1.C:449
 testUnfold1.C:450
 testUnfold1.C:451
 testUnfold1.C:452
 testUnfold1.C:453
 testUnfold1.C:454
 testUnfold1.C:455
 testUnfold1.C:456
 testUnfold1.C:457
 testUnfold1.C:458
 testUnfold1.C:459
 testUnfold1.C:460
 testUnfold1.C:461
 testUnfold1.C:462
 testUnfold1.C:463
 testUnfold1.C:464
 testUnfold1.C:465
 testUnfold1.C:466
 testUnfold1.C:467
 testUnfold1.C:468
 testUnfold1.C:469
 testUnfold1.C:470
 testUnfold1.C:471
 testUnfold1.C:472
 testUnfold1.C:473
 testUnfold1.C:474
 testUnfold1.C:475
 testUnfold1.C:476
 testUnfold1.C:477
 testUnfold1.C:478
 testUnfold1.C:479
 testUnfold1.C:480
 testUnfold1.C:481
 testUnfold1.C:482
 testUnfold1.C:483
 testUnfold1.C:484
 testUnfold1.C:485
 testUnfold1.C:486
 testUnfold1.C:487
 testUnfold1.C:488
 testUnfold1.C:489
 testUnfold1.C:490
 testUnfold1.C:491
 testUnfold1.C:492
 testUnfold1.C:493
 testUnfold1.C:494
 testUnfold1.C:495
 testUnfold1.C:496
 testUnfold1.C:497
 testUnfold1.C:498
 testUnfold1.C:499
 testUnfold1.C:500
 testUnfold1.C:501
 testUnfold1.C:502
 testUnfold1.C:503
 testUnfold1.C:504
 testUnfold1.C:505
 testUnfold1.C:506
 testUnfold1.C:507
 testUnfold1.C:508
 testUnfold1.C:509