ROOT logo
/**
 * A task to do a comparison between tracklets and clusers in the SPD
 * 
 * Since the class SPDComparisonTask derives from a compiled class
 * (AliAnalysisTaskSE) we need to compile that code.  The script will,
 * when executed in the AliROOT prompt load it self again and byte
 * compile it with the preprocessor flag BUILD defined.  \
 *
 * The preprocessor define BUILD is not defined at first, when the
 * script is loaded using 
 * 
 * @verbatim 
 *   Root> .x SPDComparison.C 
 * @endverbatim 
 * 
 * which means that CINT will only see the function SPDComparison.
 * In that function, we define the BUILD preprocessor symbol 
 *
 * @code 
 *   gSystem->AddIncludePath("-DBUILD=1 ...");
 * @endcode 
 * 
 * and then ACLic compile ourselves 
 * 
 * @code 
 *   gROOT->LoadMacro("./SPDComparison.C++");
 * @endcode 
 * 
 * But since BUILD is now defined, it means that ACLic will only see
 * the class and not the function SPDComparison. 
 * 
 * This trick hinges on that when you initially load the script and
 * when it is done inside the script it is done using two distinct
 * paths - otherwise ROOT will try to unload the script first, and
 * that fails.   
 * 
 */

#ifdef BUILD
#include <AliAnalysisTaskSE.h>
#include <AliESDEvent.h>
#include <AliLog.h>
#include <AliMultiplicity.h>
#include "AliFMDEventInspector.h"
#include "AliAODForwardMult.h"
#include <TH1D.h>
#include <TH2D.h>
#include <TGraphErrors.h>
#include <TList.h>
#include <TObjArray.h>
#include <TMath.h>


/**
 * @class SPDComparisonTask 
 *
 * Task to do forward/backward correlations using FMD and SPD
 * (cluster) data. 
 * 
 * The class contains the sub-structure Bin that represents an
 * @f$\eta@f$ range.  One can add (possibly overlapping) @f$\eta@f$
 * ranges by calling the member function AddBin 
 * 
 * @ingroup pwglf_forward_scripts
 */
class SPDComparisonTask : public AliAnalysisTaskSE
{
public:
  /**
   * A single vertex bin - contains summed histogram of clusters and
   * tracklets 
   * 
   */
  struct VtxBin : public TNamed
  {
    /** 
     * Static function to get a bin name 
     * 
     * @param low   Low limit 
     * @param high  High limit 
     * 
     * @return Pointer to static array
     */
    static const char* BinName(Double_t low, Double_t high)
    {
      static TString n;
      n = (Form("vtx%+05.1f_%+05.1f", low, high));
      n.ReplaceAll("+", "p");
      n.ReplaceAll("-", "m");
      n.ReplaceAll(".", "d");
      return n.Data();
    }
    /** 
     * Constructor
     * 
     * @param low   Low limit 
     * @param high  High limit 
     */
    VtxBin(Double_t low, Double_t high)
      : TNamed(BinName(low,high), ""),
	fClusters(0), 
	fTracklets(0)
    {
    }
    /** 
     * Constructor
     */
    VtxBin() : TNamed(), fClusters(0), fTracklets(0) {}
    /** 
     * Copy constructor 
     * 
     * @param o Object to copy from 
     */
    VtxBin(const VtxBin& o) : 
      TNamed(o), fClusters(o.fClusters), fTracklets(o.fTracklets) 
    {}
    /** 
     * Assignment operator 
     * 
     * @param o Object to assign from 
     * 
     * @return Reference to this.
     */
    VtxBin& operator=(const VtxBin& o) 
    {
      TNamed::operator=(o);
      fClusters = o.fClusters;
      fTracklets = o.fTracklets;
      return *this;
    }
    /** 
     * Define outputs
     * 
     * @param l    Output list
     * @param eta  Eta axis to use 
     */
    void DefineOutput(TList* l, const TAxis& eta)
    {
      TList* ll = new TList;
      ll->SetName(GetName());
      ll->SetOwner();
      l->Add(ll);

      fClusters = new TH1D("clusters", "Clusters", 
			   eta.GetNbins(), 
			   eta.GetXmin(), 
			   eta.GetXmax());
      fClusters->SetXTitle("#eta");
      fClusters->SetYTitle("dN/d#eta");
      fClusters->SetDirectory(0);
      fClusters->Sumw2();
      ll->Add(fClusters);

      fTracklets = static_cast<TH1D*>(fClusters->Clone("tracklets"));
      fTracklets->SetTitle("Tracklets");
      fTracklets->SetDirectory(0);
      ll->Add(fTracklets);
    }
    /** 
     * Process the input 
     * 
     * @param spdmult  Input 
     */
    void Process(const AliMultiplicity* spdmult)
    {
      //Filling clusters in layer 1 used for tracklets...
      for(Int_t j = 0; j< spdmult->GetNumberOfTracklets();j++) {
	Double_t eta = spdmult->GetEta(j);
	fClusters->Fill(eta);
	fTracklets->Fill(eta);
      }

      //...and then the unused ones in layer 1 
      for(Int_t j = 0; j< spdmult->GetNumberOfSingleClusters();j++) {
	Double_t eta = -TMath::Log(TMath::Tan(spdmult->GetThetaSingle(j)/2.));
	fClusters->Fill(eta);
      }
    }
    /** 
     * End of job processing
     * 
     * @param l        Output list
     * @param nEvents  Number of events to normalise to 
     */
    void Finish(TList* l, Int_t nEvents)
    {
      TList* ll  = static_cast<TList*>(l->FindObject(GetName()));
      fClusters  = static_cast<TH1D*>(ll->FindObject("clusters"));
      fTracklets = static_cast<TH1D*>(ll->FindObject("tracklets"));
      fClusters->Scale(1. / nEvents, "width");
      fTracklets->Scale(1. / nEvents, "width");
      TH1* ratio = static_cast<TH1D*>(fClusters->Clone("ratio"));
      ratio->SetDirectory(0);
      ratio->SetTitle("Clusters/Tracklets");
      ratio->Divide(fTracklets);
      ll->Add(ratio);
    }
    TH1D* fClusters;  // Cluster histogram
    TH1D* fTracklets; // Tracklet histogram 
  };
    
  //__________________________________________________________________
  /** 
   * Default constructor - do not use 
   */
  SPDComparisonTask() 
    : AliAnalysisTaskSE(), 
      fInspector(),
      fBins(0),
      fEvents(0),
      fList(0), 
      fFirstEvent(true), 
      fVertexAxis(20, -20, 20)
  {}
  /** 
   * Constructor with a single argument. 
   */
  SPDComparisonTask(const char*) 
    : AliAnalysisTaskSE("spdComparision"), 
      fInspector("eventInspector"),
      fBins(0),
      fEvents(0),
      fList(0), 
      fFirstEvent(true), 
      fVertexAxis(20, -20, 20)
  {
    // Declare our output container 
    DefineOutput(1,TList::Class());
  }
  /** 
   * Copy constructor 
   * 
   * @param o Object to copy from 
   */
  SPDComparisonTask(const SPDComparisonTask& o) 
    : AliAnalysisTaskSE(o),
      fInspector(o.fInspector),
      fBins(o.fBins),
      fEvents(o.fEvents),
      fList(o.fList), 
      fFirstEvent(o.fFirstEvent), 
      fVertexAxis(20, -20, 20)
  {
    SetVertexAxis(o.fVertexAxis);
  }
  /** 
   * Assigment operator 
   * 
   * @param o Object to assign from 
   * 
   * @return Reference to this 
   */
  SPDComparisonTask& operator=(const SPDComparisonTask& o) 
  { 
    AliAnalysisTaskSE::operator=(o);
    fInspector = o.fInspector;
    fEvents    = o.fEvents;
    fList      = o.fList;
    fFirstEvent = o.fFirstEvent;
    SetVertexAxis(o.fVertexAxis);
    return *this; 
  }
  /** 
   * Destructor 
   */
  ~SPDComparisonTask() {}
  /** 
   * Set the vertex axis to use 
   * 
   * @param a Axis to set from 
   */
  void SetVertexAxis(const TAxis& a)
  {
    SetVertexAxis(a.GetNbins(), 
		  a.GetXmin(),
		  a.GetXmax());
  }
  /** 
   * Set the vertex axis to use 
   * 
   * @param n     Number of bins
   * @param xmin  Least @f$v_z@f$
   * @param xmax  Most @f$v_z@f$
   */
  void SetVertexAxis(Int_t n, Double_t xmin, Double_t xmax)
  {
    fVertexAxis.Set(n, xmin, xmax);
  }
  /** 
   * Create output objects 
   * 
   */
  void UserCreateOutputObjects()
  {
    fList = new TList;
    fList->SetName(GetName());
    fList->SetOwner();


    fEvents = new TH1D("events", "Events", 
		       fVertexAxis.GetNbins(), 
		       fVertexAxis.GetXmin(), 
		       fVertexAxis.GetXmax());
    fEvents->SetDirectory(0);
    fEvents->SetXTitle("v_{z} [cm]");
    fList->Add(fEvents);

    TAxis& vtxAxis = fVertexAxis;
    fBins = new TObjArray(vtxAxis.GetNbins(),1);
    
    TAxis etaAxis(120, -3, 3);
    for (Int_t i = 1; i <= vtxAxis.GetNbins(); i++) {
      VtxBin* bin = new VtxBin(vtxAxis.GetBinLowEdge(i), 
			       vtxAxis.GetBinUpEdge(i));
      bin->DefineOutput(fList, etaAxis);
      fBins->AddAt(bin,i);
    }


    fInspector.DefineOutput(fList);
    fInspector.Init(*(fEvents->GetXaxis()));
    
    PostData(1, fList);
  }
  /** 
   * Process one event
   * 
   * @param option Not used 
   */    
  void UserExec(Option_t* option="")
  {
    // Get the input data - ESD event
    AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
    if (!esd) { 
      AliWarning("No ESD event found for input event");
      return;
    }
    
    //--- Read run information -----------------------------------------
    if (fFirstEvent && esd->GetESDRun()) {
      fInspector.ReadRunDetails(esd);
      
      AliInfo(Form("Initializing with parameters from the ESD:\n"
		   "         AliESDEvent::GetBeamEnergy()   ->%f\n"
		   "         AliESDEvent::GetBeamType()     ->%s\n"
		   "         AliESDEvent::GetCurrentL3()    ->%f\n"
		   "         AliESDEvent::GetMagneticField()->%f\n"
		   "         AliESDEvent::GetRunNumber()    ->%d\n",
		   esd->GetBeamEnergy(),
		   esd->GetBeamType(),
		   esd->GetCurrentL3(),
		   esd->GetMagneticField(),
		   esd->GetRunNumber()));
      
      Print();
      fFirstEvent = false;
    }
    
    // Some variables 
    UInt_t   triggers; // Trigger bits
    Bool_t   lowFlux;  // Low flux flag
    UShort_t iVz;      // Vertex bin from ESD
    Double_t vZ;       // Z coordinate from ESD
    Double_t cent;     // Centrality 
    UShort_t nClusters;// Number of SPD clusters 
    // Process the data 
    UInt_t retESD = fInspector.Process(esd, triggers, lowFlux, iVz, vZ, 
				       cent, nClusters);
    Bool_t isInel   = triggers & AliAODForwardMult::kB;
    Bool_t hasVtx   = retESD == AliFMDEventInspector::kOk;

    if (!isInel || !hasVtx) return;
    fEvents->Fill(vZ);

    const AliMultiplicity* spdmult = esd->GetMultiplicity();

    VtxBin* bin = static_cast<VtxBin*>(fBins->At(iVz));
    if (!bin) { 
      AliError(Form("No bin @ %d (%fcm)", iVz, vZ));
      return;
    }
    bin->Process(spdmult);
    
    // Post data to output container 
    PostData(1,fList);
  }
  /** 
   * Called at the end of the processing 
   * 
   * @param option 
   */
  void Terminate(Option_t* option="")
  {
    fList = dynamic_cast<TList*>(GetOutputData(1));
    if (!fList) {
      AliError("No output list defined");
      return;
    }

    fEvents    = static_cast<TH1D*>(fList->FindObject("events"));

    Int_t nEvents = fEvents->GetEntries();
    AliInfo(Form("Got a total of %d events", nEvents));

    TH1* clusters = 0;
    TH1* tracklets = 0;
    TIter next(fBins);
    VtxBin* bin = 0;
    Int_t i = 1;
    while ((bin = static_cast<VtxBin*>(next()))) { 
      bin->Finish(fList, fEvents->GetBinContent(i++));
      if (!clusters) clusters = static_cast<TH1D*>(bin->fClusters->Clone());
      else clusters->Add(bin->fClusters);
      if (!tracklets) tracklets = static_cast<TH1D*>(bin->fTracklets->Clone());
      else tracklets->Add(bin->fTracklets);
    }
    clusters->SetDirectory(0);
    tracklets->SetDirectory(0);
    clusters->Scale(1. / i);
    tracklets->Scale(1. / i);
    

    TH1D* ratio = static_cast<TH1D*>(clusters->Clone("ratio"));
    ratio->SetDirectory(0);
    ratio->SetTitle("Clusters/Tracklets");
    ratio->Divide(tracklets);

    fList->Add(clusters);
    fList->Add(tracklets);
    fList->Add(ratio);
  }
protected: 
  AliFMDEventInspector fInspector; // Inspector
  TObjArray*           fBins;      // Vertex bins 
  TH1D*                fEvents;    // Events 
  TList*               fList;      // Output list
  Bool_t               fFirstEvent;// First event flag 
  TAxis                fVertexAxis;

  ClassDef(SPDComparisonTask,1); // BF analysis
};
#else

//====================================================================
/** 
 * Run the analysis 
 * 
 * @param esddir  Input directory 
 * @param nEvents Number of events, negative means all 
 */
void
SPDComparison(const char* esddir, Int_t nEvents=-1)
{
  // --- Libraries to load -------------------------------------------
  gROOT->Macro("$ALICE_ROOT/PWGLF/FORWARD/analysis2/scripts/LoadLibs.C");

  // --- Our data chain ----------------------------------------------
  gROOT->LoadMacro("$ALICE_ROOT/PWGLF/FORWARD/analysis2/scripts/MakeChain.C");
  TChain* chain = MakeChain("ESD", esddir, true);
  // If 0 or less events is select, choose all 
  if (nEvents <= 0) nEvents = chain->GetEntries();
  
  // --- Manager -----------------------------------------------------
  AliAnalysisManager* mgr = new AliAnalysisManager("FB", "FB train");
  AliAnalysisManager::SetCommonFileName("spd_comps.root");

  // --- AOD input handler -------------------------------------------
  AliESDInputHandler *inputHandler = new AliESDInputHandler();
  mgr->SetInputEventHandler(inputHandler);      
   
  // --- compile our code --------------------------------------------
  gSystem->AddIncludePath("-I${ALICE_ROOT}/PWGLF/FORWARD/analysis2 "
                          "-I${ALICE_ROOT}/ANALYSIS "
                          "-I${ALICE_ROOT}/include -DBUILD=1");
  gROOT->LoadMacro("./SPDComparison.C++g");
  
  // --- Make our object ---------------------------------------------
  SPDComparisonTask* task = new SPDComparisonTask("SPD_COMP");
  task->SetVertexAxis(10, -10, 10);
  mgr->AddTask(task);

  // --- create containers for input/output --------------------------
  AliAnalysisDataContainer *sums = 
    mgr->CreateContainer("spdComp", TList::Class(), 
                         AliAnalysisManager::kOutputContainer, 
                         AliAnalysisManager::GetCommonFileName());
  // --- connect input/output ----------------------------------------
  mgr->ConnectInput(task, 0, mgr->GetCommonInputContainer());
  mgr->ConnectOutput(task, 1, sums);

  // --- Run the analysis --------------------------------------------
  TStopwatch t;
  if (!mgr->InitAnalysis()) {
    Error("SPDComparison", "Failed to initialize analysis train!");
    return;
  }
  // Some informative output 
  mgr->PrintStatus();
  // mgr->SetDebugLevel(3);
  if (mgr->GetDebugLevel() < 1) mgr->SetUseProgressBar(kTRUE,100);

  // Run the train 
  t.Start();
  Printf("=== RUNNING ANALYSIS on %9d events ==================", nEvents);
  mgr->StartAnalysis("local", chain, nEvents);
  t.Stop();
  t.Print();

  DrawSPDComparison();
}

//====================================================================
/** 
 * Draw results
 * 
 * @param filename 
 */
void
DrawSPDComparison(const char* filename="spd_comps.root")
{
  // --- Open the file -----------------------------------------------
  TFile* file = TFile::Open(filename, "READ");
  if (!file) { 
    Error("DrawSPDComparison", "Failed to open file %s", filename);
    return;
  }
  
  // --- Get our list ------------------------------------------------
  TList* spd = static_cast<TList*>(file->Get("spdComp"));
  if (!spd) { 
    Error("DrawSPDComparison", "Failed to get list SPD_COMP from file %s",
	  filename);
    return;
  }

  // --- Loop over list and get correlation plots --------------------
  TH1* clusters  = static_cast<TH1*>(spd->FindObject("clusters"));
  TH1* tracklets = static_cast<TH1*>(spd->FindObject("tracklets"));
  TH1* ratio     = static_cast<TH1*>(spd->FindObject("ratio"));
  TH1* events    = static_cast<TH1*>(spd->FindObject("events"));

  // --- Set style parameters ----------------------------------------
  gStyle->SetPalette(1);
  gStyle->SetTitleFillColor(0);
  gStyle->SetTitleStyle(0);
  gStyle->SetTitleBorderSize(0);
  gStyle->SetOptStat(0);
  gStyle->SetTitleX(0.69);
  gStyle->SetTitleY(0.99);
  gStyle->SetTitleW(0.30);
  gStyle->SetTitleH(0.10);

  // --- Make canvas for correlation plots and plot them -------------
  TCanvas* c1 = new TCanvas("c1", "c1", 900, 600);
  c1->SetTopMargin(0.05);
  c1->SetRightMargin(0.05);
  c1->SetFillColor(0);
  c1->SetFillStyle(0);
  c1->SetBorderSize(0);
  c1->SetBorderMode(0);
  
  c1->cd();
  TPad* p1 = new TPad("p1","p1", 0.0, 0.3, 1.0, 1.0, 0, 0, 0);
  p1->SetFillColor(0);
  p1->SetFillStyle(0);
  p1->SetBorderSize(0);
  p1->SetBorderMode(0);
  p1->SetTopMargin(0.05);
  p1->SetBottomMargin(0.001);
  p1->SetRightMargin(0.05);
  p1->Draw();
  p1->cd();
  clusters->SetMarkerStyle(20);
  clusters->SetMarkerColor(kRed+1);
  clusters->Draw();

  tracklets->SetMarkerStyle(20);
  tracklets->SetMarkerColor(kBlue+1);
  tracklets->Draw("same");

  TString v(Form("%+5.1f<|v_{z}|<%+5.1f",
		 events->GetXaxis()->GetXmin(),
		 events->GetXaxis()->GetXmax()));
  TLatex* ltx = new TLatex(.2,.80, v.Data());
  ltx->Draw();

  TLegend* l = p1->BuildLegend(.6,.6,.94,.94);
  l->SetFillColor(0);
  l->SetFillStyle(0);
  l->SetBorderSize(0);

  c1->cd();
  TPad* p2 = new TPad("p2","p2", 0.0, 0.0, 1.0, 0.3, 0, 0, 0);
  p2->SetFillColor(0);
  p2->SetFillStyle(0);
  p2->SetBorderSize(0);
  p2->SetBorderMode(0);
  p2->SetTopMargin(0.001);
  p2->SetBottomMargin(0.2);
  p2->SetRightMargin(0.05);
  p2->Draw();
  p2->cd();
  ratio->SetMarkerStyle(20);
  ratio->SetMarkerColor(kGray+1);
  ratio->Draw();
  ratio->GetYaxis()->SetRangeUser(0,2);

  c1->SaveAs("SPDComparison.png");
  c1->SaveAs("SPDComparison.root");
}

    
  
  
  
#endif
//
// EOF
//

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