Ignore:
Timestamp:
01/20/05 10:09:24 (20 years ago)
Author:
domingo
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars/mtemp/mifae/library
Files:
2 added
2 deleted
8 edited

Legend:

Unmodified
Added
Removed

  • trunk/MagicSoft/Mars/mtemp/mifae/library/IFAELinkDef.h

    r5672 r5904  
    2727#pragma link C++ class MTopologyCalc+;
    2828#pragma link C++ class MImageParDisp+;
    29 #pragma link C++ class MDisp+;
     29#pragma link C++ class MDispParameters+;
    3030#pragma link C++ class MDispCalc+;
    3131#pragma link C++ class MHDisp+;

  • trunk/MagicSoft/Mars/mtemp/mifae/library/MDispCalc.cc

    r5671 r5904  
    2828//   MDispCalc                                                             //
    2929//                                                                         //
    30 //   This task calculates Disp                                             //
    31 //   (parameters are taken from the container MDisp)                       //
    32 //   Also the matrix of variables to be used in the Disp optimization      //
    33 //   is defined                                                            //
     30//   This task calculates Disp with a given parameterization               //
     31//   (parameters are stored in the MDispParameters container)              //
     32//   Also the matrix of variables to be used in the Disp                   //
     33//   parameterization is defined                                           //
    3434//                                                                         //
    3535//                                                                         //
    … …  
    4848#include "MPointingPos.h"
    4949#include "MImageParDisp.h"
    50 #include "MDisp.h"
     50#include "MDispParameters.h"
    5151
    5252#include "MHMatrix.h"
    … …  
    6262using namespace std;
    6363
     64static const Int_t nPars=5;   // number of parameters used in the Disp expression
     65
     66enum dispVar_t {kSize,kWidth,kLength,kConc,kLeakage1,kLeakage2,kTotVar};  // enum variables for the
     67                                                                          // matrix columns mapping
     68
    6469// --------------------------------------------------------------------------
    6570//
    6671// constructor
    6772//
    68 MDispCalc::MDispCalc(const char *imagepardispname, const char *dispname,
     73MDispCalc::MDispCalc(const char *imagepardispname,
     74                     const char *dispparametersname,
    6975                     const char *name, const char *title)
    70   : fImageParDispName(imagepardispname), fDispName(dispname)
     76  :     fImageParDispName(imagepardispname),
     77        fDispParametersName(dispparametersname),
     78        fLogSize0(3.), fLength0(0.1), fWidth0(0.05),
     79        fConc0(0.35), fLeakage10(0.05), fLeakage20(0.1)
    7180{
    7281    fName  = name  ? name  : "MDispCalc";
    … …  
    7685
    7786    // initialize mapping array dimension to the number of columns of fMatrix
    78     fMap.Set(18);
     87    fMap.Set(kTotVar);
    7988}
    8089
    … …  
    106115    }
    107116
    108     fDisp = (MDisp*)pList->FindObject(fDispName, "MDisp");
    109     if (!fDisp)
    110     {
    111         *fLog << err << fDispName << " [MDisp] not found... aborting." << endl;
    112         return kFALSE;
    113     }
     117    fDispParameters = (MDispParameters*)pList->FindObject(fDispParametersName, "MDispParameters");
     118    if (!fDispParameters)
     119    {
     120        *fLog << err << fDispParametersName << " [MDispParameters] not found... aborting." << endl;
     121        return kFALSE;
     122    }
     123    fDispParameters->GetParameters().Set(nPars);
     124    fDispParameters->GetStepsizes().Set(nPars);
    114125
    115126
    … …  
    160171    if (!fPointing)
    161172    {
    162         *fLog << err << "MPointingPos not found... "
    163               << "MMcEvt is going to be used to get Theta and Phi."
    164               << endl;
    165         //        return kFALSE;
     173        *fLog << err << "MPointingPos not found... aborting." << endl;
     174        return kFALSE;
    166175    }
    167176
    … …  
    175184//
    176185// You can use this function if you want to use a MHMatrix instead of the
    177 // given containers for the Disp optimization. This function adds all
     186// given containers for the Disp parameterization. This function adds all
    178187// necessary columns to the given matrix. Afterwards, you should fill
    179188// the matrix with the corresponding data (eg from a file by using
    180189// MHMatrix::Fill). Then, if you loop through the matrix (eg using
    181 // MMatrixLoop), MDispCalc::Process will take the values from the matrix
     190// MMatrixLoop), MDispCalc::Calc will take the values from the matrix
    182191// instead of the containers.
    183192//
    … …  
    189198    fMatrix = mat;
    190199
    191     fMap[0]  = fMatrix->AddColumn("MSrcPosCam.fX");
    192     fMap[1]  = fMatrix->AddColumn("MSrcPosCam.fY");
    193     fMap[2]  = fMatrix->AddColumn("MHillas.fMeanX");
    194     fMap[3]  = fMatrix->AddColumn("MHillas.fMeanY");
    195     fMap[4]  = fMatrix->AddColumn("MHillas.fDelta");
    196 
    197     fMap[5]  = fMatrix->AddColumn("MHillas.fSize");
    198     fMap[6]  = fMatrix->AddColumn("MHillas.fWidth");
    199     fMap[7]  = fMatrix->AddColumn("MHillas.fLength");
    200 
    201     fMap[8]  = fMatrix->AddColumn("MNewImagePar.fConc");
    202     fMap[9]  = fMatrix->AddColumn("MNewImagePar.fLeakage1");
    203     fMap[10] = fMatrix->AddColumn("MNewImagePar.fLeakage2");
    204 
    205     fMap[11] = fMatrix->AddColumn("MHillasExt.fM3Long");
    206     fMap[12] = fMatrix->AddColumn("MHillasExt.fAsym");
    207 
    208     fMap[13]  = fMatrix->AddColumn("MMcEvt.fEnergy");
    209     fMap[14]  = fMatrix->AddColumn("MMcEvt.fImpact");
    210     fMap[15]  = fMatrix->AddColumn("MMcEvt.fLongitmax");
    211    
    212     if (fPointing)
    213     {
    214       fMap[16]  = fMatrix->AddColumn("MPointingPos.fZd");
    215       fMap[17]  = fMatrix->AddColumn("MPointingPos.fAz");
    216     }
    217     else
    218     {
    219       fMap[16]  = fMatrix->AddColumn("MMcEvt.fTelescopeTheta*kRad2Deg");
    220       fMap[17]  = fMatrix->AddColumn("MMcEvt.fTelescopePhi*kRad2Deg");
    221     }
    222 }
    223 
    224 
    225 // --------------------------------------------------------------------------
    226 //
    227 // Returns the Matrix and the mapped value for each Matrix column
    228 //
    229 MHMatrix* MDispCalc::GetMatrixMap(TArrayI &map)
    230 {
    231     map = fMap;
    232 
    233     return fMatrix;
     200    fMap[kSize]      = fMatrix->AddColumn("MHillas.fSize");
     201    fMap[kWidth]     = fMatrix->AddColumn("MHillas.fWidth");
     202    fMap[kLength]    = fMatrix->AddColumn("MHillas.fLength");
     203
     204    fMap[kConc]      = fMatrix->AddColumn("MNewImagePar.fConc");
     205    fMap[kLeakage1]  = fMatrix->AddColumn("MNewImagePar.fLeakage1");
     206    fMap[kLeakage2]  = fMatrix->AddColumn("MNewImagePar.fLeakage2");
    234207}
    235208
    … …  
    258231{
    259232    // get variables needed to compute disp from the matrix or the containers
    260     const Double_t size     = fMatrix ? GetVal(5)  : fHil->GetSize();
    261     const Double_t width0   = fMatrix ? GetVal(6)  : fHil->GetWidth();
    262     const Double_t length0  = fMatrix ? GetVal(7)  : fHil->GetLength();
    263     const Double_t conc     = fMatrix ? GetVal(8)  : fNewPar->GetConc();
    264     const Double_t leakage1 = fMatrix ? GetVal(9) : fNewPar->GetLeakage1();
    265     const Double_t leakage2 = fMatrix ? GetVal(10) : fNewPar->GetLeakage2();
     233    const Double_t size     = fMatrix ? GetVal(kSize)     : fHil->GetSize();
     234    const Double_t width0   = fMatrix ? GetVal(kWidth)    : fHil->GetWidth();
     235    const Double_t length0  = fMatrix ? GetVal(kLength)   : fHil->GetLength();
     236    const Double_t conc     = fMatrix ? GetVal(kConc)     : fNewPar->GetConc();
     237    const Double_t leakage1 = fMatrix ? GetVal(kLeakage1) : fNewPar->GetLeakage1();
     238    const Double_t leakage2 = fMatrix ? GetVal(kLeakage2) : fNewPar->GetLeakage2();
    266239
    267240    // convert to deg
    … …  
    269242    const Double_t length  = length0 * fMm2Deg;
    270243
    271     // create an array to pass the variables to MDisp::Calc
    272     Int_t numimagevars = 6; 
    273     TArrayD imagevars(numimagevars);
    274     imagevars[0]  = log10(size);
    275     imagevars[1]  = width;
    276     imagevars[2]  = length;
    277     imagevars[3]  = conc;
    278     imagevars[4]  = leakage1;
    279     imagevars[5]  = leakage2;
     244    // create an array to pass the variables to MDispCalc::Calc
     245    TArrayD imagevars(kTotVar);
     246    imagevars[kSize]      = log10(size);
     247    imagevars[kWidth]     = width;
     248    imagevars[kLength]    = length;
     249    imagevars[kConc]      = conc;
     250    imagevars[kLeakage1]  = leakage1;
     251    imagevars[kLeakage2]  = leakage2;
    280252
    281253    // compute Disp
    282     Double_t disp = fDisp->Calc(imagevars);
     254    Double_t disp = Calc(imagevars);
    283255
    284256    // save value into MImageParDisp container
    … …  
    288260    return kTRUE;
    289261}
     262
     263
     264// --------------------------------------------------------------------------
     265//
     266// Set the Disp parameterization and Calculate Disp
     267//
     268//
     269Double_t MDispCalc::Calc(TArrayD &imagevar)
     270{
     271    // get parameters
     272    const TArrayD& p = fDispParameters->GetParameters();
     273
     274    // get image variables to be used in the Disp parameterization
     275    Double_t logsize  = imagevar[0];
     276    Double_t width    = imagevar[1];
     277    Double_t length   = imagevar[2];
     278    //    Double_t conc     = imagevar[3];
     279    //    Double_t leakage1 = imagevar[4];
     280    //    Double_t leakage2 = imagevar[5];
     281   
     282    // Disp parameterization to be optimized
     283    //  Note: fLogSize0, fLength0... variables are introduced to uncorrelate as much
     284    //        as possible the parameters in the Disp expression, with the purpose of
     285    //        helping the minimization algorithm to converge. They are set approx.
     286    //        to their distribution mean value in the MDisp constructor, but can be
     287    //        changed using the corresponding set function.
     288    //
     289
     290    //    Double_t disp = p[0] + p[1]*(logsize-fLogSize0)
     291    //      + (p[2] + p[3]*(logsize-fLogSize0))*width/length;
     292
     293    Double_t disp = p[0] + p[1]*(logsize-fLogSize0) + p[4]*(length-fLength0)
     294      + (p[2] + p[3]*(logsize-fLogSize0))*width/length;
     295
     296    //    Double_t disp = p[0] + p[1]*(logsize-fLogSize0) + p[4]*(length-fLength0)
     297    //      + (p[2] + p[3]*(logsize-fLogSize0))*length/width;
     298
     299    //    Double_t disp = p[0] + p[1]*(logsize-fLogSize0) + p[4]*(length-fLength0)
     300    //      + (p[2] + p[3]*(logsize-fLogSize0) + p[5]*(length-fLength0))*width/length;
     301
     302    //    Double_t disp = p[0] + p[1]*(logsize-fLogSize0) + p[4]*(width-fWidth0)
     303    //      + (p[2] + p[3]*(logsize-fLogSize0))*width/length;   // + p[5]*(width-fWidth0))*width/length;
     304
     305    //    Double_t disp = p[0] + p[1]*(logsize-fLogSize0) + p[4]*(conc-fConc0)
     306    //      + (p[2] + p[3]*(logsize-fLogSize0))*width/length;   // + p[5]*(conc-fConc0))*width/length;
     307
     308    //    Double_t disp = ( p[0] + p[1]*(logsize-fLogSize0)
     309    //                + p[2]*pow(logsize-fLogSize0,2)
     310    //                + p[3]*pow(logsize-fLogSize0,3)
     311    //                + p[4]*pow(logsize-fLogSize0,4) )
     312    //                    *( 1/(1+width/length) );
     313
     314    //    Double_t disp = ( p[0] + p[1]*(logsize-fLogSize0)
     315    //                + p[2]*pow(logsize-fLogSize0,2)
     316    //                + p[3]*pow(logsize-fLogSize0,3)
     317    //                + p[4]*pow(logsize-fLogSize0,4) )
     318    //                    *(1-width/length);
     319
     320    return disp;
     321}
     322
    290323//===========================================================================
    291324

  • trunk/MagicSoft/Mars/mtemp/mifae/library/MDispCalc.h

    r5671 r5904  
    2121class MPointingPos;
    2222class MImageParDisp;
    23 class MDisp;
     23class MDispParameters;
    2424class MHMatrix;
    2525class MParList;
    … …  
    2929private:
    3030
    31     MSrcPosCam    *fSrcPos;
    32     MHillas       *fHil;
    33     MHillasExt    *fHilExt;
    34     MNewImagePar  *fNewPar;
    35     MMcEvt        *fMcEvt;
    36     MPointingPos  *fPointing;      // NOTE: take (theta,phi) from this container
    37                                    // when it is available (MC calibrated files
    38                                    // or data files)
    39 
    40     MImageParDisp *fImageParDisp;  //! output container for Disp
    41     MDisp         *fDisp;          //  container for Disp parameters
     31    MSrcPosCam      *fSrcPos;
     32    MHillas         *fHil;
     33    MHillasExt      *fHilExt;
     34    MNewImagePar    *fNewPar;
     35    MMcEvt          *fMcEvt;
     36    MPointingPos    *fPointing;       
     37   
     38    MImageParDisp   *fImageParDisp;    //! output container for Disp value
     39    MDispParameters *fDispParameters;  //!  container for Disp parameters
    4240
    4341    TString  fImageParDispName;    // name of container to store Disp
    44     TString  fDispName;            // name of container for Disp parameters
     42    TString  fDispParametersName;  // name of container for Disp parameters
     43
     44    Double_t fLogSize0;            // Variables introduced in the Disp expression
     45    Double_t fLength0;             // to shift the minimization around the mean
     46    Double_t fWidth0;              // values of the variables, so it makes easier
     47    Double_t fConc0;               // to MINUIT to converge. Default values are set
     48    Double_t fLeakage10;           // in the constructor (to their standard mean
     49    Double_t fLeakage20;           // distribution values) but can be changed with
     50                                   // the corresponding set function.
    4551
    4652    Double_t fMm2Deg;              // conversion factor from mm to deg
    4753
    4854    MHMatrix *fMatrix;             // matrix defined to have as much columns as
    49                                    // variables wanted to be accessible during
    50                                    // Disp optimization, and as much rows as
    51                                    // events used for the optimization
    52                                    // (filled at MFindDisp)
     55                                   // variables wanted to parameterize Disp,
     56                                   // and as much rows as events used for the
     57                                   // optimization (filled at MFindDisp)
    5358
    54     TArrayI  fMap;                 // array to store the matrix mapping column
    55                                    // numbers corresponding to each added variable
     59    TArrayI  fMap;                 // array storing the matrix mapping column
     60                                   // numbers corresponding to each variable
     61                                   // added to fMatrix
    5662
    5763
    … …  
    5965    Int_t Process();
    6066
    61     Double_t GetVal(Int_t i) const; // get value of variable entered in the matrix
    62                                     // at column fMap[i]
     67    Double_t GetVal(Int_t i) const;        // get value of variable entered in the
     68                                           // matrix at column fMap[i]
    6369
    6470
    … …  
    7076
    7177    void InitMapping(MHMatrix *mat);       // define the matrix of variables
    72     MHMatrix* GetMatrixMap(TArrayI &map);  // get matrix and its mapping array
     78                                           // needed for the Disp parameterization
     79
     80    Double_t Calc(TArrayD &imagevar);      // calculate Disp with a given expression
     81
     82    void SetLogSize0(Double_t newmeanval)  { fLogSize0  = newmeanval; }
     83    void SetWidth0(Double_t newmeanval)    { fWidth0    = newmeanval; }
     84    void SetLength0(Double_t newmeanval)   { fLength0   = newmeanval; }
     85    void SetConc0(Double_t newmeanval)     { fConc0     = newmeanval; }
     86    void SetLeakage10(Double_t newmeanval) { fLeakage10 = newmeanval; }
     87    void SetLeakage20(Double_t newmeanval) { fLeakage20 = newmeanval; }
    7388
    7489    ClassDef(MDispCalc, 0) // Task to evaluate Disp

  • trunk/MagicSoft/Mars/mtemp/mifae/library/MFindDisp.cc

    r5671 r5904  
    4848#include "MBinning.h"
    4949#include "MContinue.h"
    50 #include "MDisp.h"
     50#include "MDispParameters.h"
    5151#include "MDispCalc.h"
    5252#include "MDataElement.h"
    … …  
    7272#include "MMatrixLoop.h"
    7373#include "MMinuitInterface.h"
     74#include "MParameters.h"
    7475#include "MParList.h"
    7576#include "MProgressBar.h"
    … …  
    9596//
    9697// - the parameters to be varied in the minimization are the parameters
    97 //   appearing in the parametrization of Disp (MDisp::Calc())
     98//   appearing in the parametrization of Disp (MDispCalc::Calc())
    9899//
    99100//------------------------------------------------------------------------
    … …  
    104105    cout <<  "entry fcnDisp" << endl;
    105106
    106     // save pointer to the MINUIT oject (for optimizing Disp) for the case
     107    // save pointer to the MINUIT object (for optimizing Disp) for the case
    107108    // that it is overwritten by the pointer to another Minuit object
    108109    //    TMinuit *savePointer = gMinuit;
    … …  
    116117
    117118    // get needed Disp containers from the existing parList
    118     MDisp *disp = (MDisp*)plistfcn->FindObject("MDisp");
    119     if (!disp)
    120     {
    121         gLog << "fcnDisp : MDisp object '" << "MDisp"
     119    MDispParameters *dispparams = (MDispParameters*)plistfcn->FindObject("MDispParameters");
     120    if (!dispparams)
     121    {
     122        gLog << "fcnDisp : MDispParameters object '" << "MDispParameters"
    122123            << "' not found... aborting" << endl;
    123124        return;
    … …  
    132133    }
    133134
     135    MParameterD *minpar = (MParameterD*)plistfcn->FindObject("MinimizationParameter");
     136
    134137    //
    135     // transfer current parameter values to MDisp
     138    // transfer current parameter values to MDispParameters
    136139    //
    137140    // Attention : npar is the number of variable parameters
    … …  
    142145    gMinuit->mnstat(fMin, fEdm, fErrdef, fNpari, fNparx, fIstat);
    143146
    144     disp->SetParameters(TArrayD(fNparx, par));
     147    dispparams->SetParameters(TArrayD(fNparx, par));
    145148
    146149    // checking that parameters have been properly set
    147     TArrayD checkparameters = disp->GetParameters();
     150    TArrayD checkparameters = dispparams->GetParameters();
    148151    gLog << "fcnDisp : fNpari, fNparx =" << fNpari << ",  "
    149152         << fNparx  << endl;
    … …  
    155158
    156159    //
    157     // loop over the events in the training matrix to compute the Chi^2
    158     // for the current parameter values
     160    // loop over the events in the training matrices to compute the parameter
     161    // to minimize for the current Disp parameter values
    159162    evtloopfcn->Eventloop();
    160163
    … …  
    162165
    163166    //-------------------------------------------
    164     // get the Chi^2 value for the current values of the parameters
    165     f = hdisp->GetChi2();
     167    // get the Minimization parameter value for the current values
     168    // of the Disp parameters
     169    if (minpar)
     170      f = minpar->GetVal();
     171    else
     172    {
     173      gLog << "fcnDisp : MParameterD object '" << "MinimizationParameter"
     174           << "' not found... aborting" << endl;
     175      return;
     176    }
    166177}
    167178
    … …  
    184195    // (done at MFDisp constructor)
    185196    fDispFilter = fdisp;
    186     *fLog << inf << "MFindDisp::MFindDisp; address of MFDisp object = "
    187           << fDispFilter << endl;   
    188197
    189198    // matrices to contain the training/test samples
    190     fMatrixTrain = new MHMatrix("MatrixTrain");
    191     fMatrixTest  = new MHMatrix("MatrixTest");
    192 
    193     // objects of MDispCalc to which these matrices are attached
     199    fMatrixTrainCalc  = new MHMatrix("MatrixTrainCalc");
     200    fMatrixTrainHists = new MHMatrix("MatrixTrainHists");
     201    fMatrixTestCalc   = new MHMatrix("MatrixTestCalc");
     202    fMatrixTestHists  = new MHMatrix("MatrixTestHists");
     203
     204    // objects of MDispCalc where the first part of the matrices mapping is defined
    194205    fDispCalcTrain = new MDispCalc("DispTrain");
    195206    fDispCalcTest  = new MDispCalc("DispTest");
    196207
     208    // objects of MHDisp where the second part of the matrices mapping is defined
     209    fHDispTrain = new MHDisp("DispTrain");
     210    fHDispTest  = new MHDisp("DispTest");
     211
    197212    // define columns of matrices
    198     fDispCalcTrain->InitMapping(fMatrixTrain);
    199     fDispCalcTest->InitMapping(fMatrixTest);
     213    // Train matrix
     214    fDispCalcTrain->InitMapping(fMatrixTrainCalc);
     215    fHDispTrain->InitMapping(fMatrixTrainHists);
     216    // Test matrix
     217    fDispCalcTest->InitMapping(fMatrixTestCalc);
     218    fHDispTest->InitMapping(fMatrixTestHists);
    200219}
     220
    201221
    202222// --------------------------------------------------------------------------
    … …  
    207227{
    208228    delete fCam;
    209     delete fMatrixTrain;
    210     delete fMatrixTest;
     229    delete fMatrixTrainCalc;
     230    delete fMatrixTrainHists;
     231    delete fMatrixTestCalc;
     232    delete fMatrixTestHists;
    211233    delete fDispCalcTrain;
    212234    delete fDispCalcTest;
     235    delete fHDispTrain;
     236    delete fHDispTest;
    213237}
    214238
    … …  
    216240// --------------------------------------------------------------------------
    217241//
    218 // Define the matrix 'fMatrixTrain' for the TRAINING sample
     242// Define the matrices 'fMatrixTrainCalc' and 'fMatrixTrainHists'
     243// for the TRAINING sample
    219244//
    220245// alltogether 'howmanytrain' events are read from file 'nametrain';
    … …  
    223248//
    224249Bool_t MFindDisp::DefineTrainMatrix(
    225                           const TString &nametrain, MH3 &hreftrain,
    226                           const Int_t howmanytrain, const TString &filetrain,
     250                          const TString &nametrain, MH3 &hreftrain, 
     251                          const Int_t howmanytrain, const TString &filetrain,
    227252                          Int_t iflag)
    228253{
    … …  
    231256
    232257    *fLog   << "=============================================" << endl;
    233     *fLog   << "Fill TRAINING Matrix from file '" << nametrain << endl;
     258    *fLog   << "Fill TRAINING Matrices from file '" << nametrain << endl;
    234259    *fLog   << "     select " << howmanytrain << " events " << endl;
    235260    if (!hreftrain.GetHist().GetEntries()==0)
    … …  
    266291    seltrain.SetName("selectTrain");
    267292
    268     MFillH filltrain(fMatrixTrain);
    269     filltrain.SetFilter(&seltrain);
    270     filltrain.SetName("fillMatrixTrain");
     293    MFillH filltraincalc(fMatrixTrainCalc);
     294    filltraincalc.SetFilter(&seltrain);
     295    filltraincalc.SetName("fillMatrixTrainCalc");
     296
     297    MFillH filltrainhists(fMatrixTrainHists);
     298    filltrainhists.SetFilter(&seltrain);
     299    filltrainhists.SetName("fillMatrixTrainHists");
    271300
    272301    //******************************
    … …  
    275304    plist.AddToList(&tlist);
    276305    plist.AddToList(fCam);
    277     plist.AddToList(fMatrixTrain);
     306    plist.AddToList(fMatrixTrainCalc);
     307    plist.AddToList(fMatrixTrainHists);
    278308
    279309    //******************************
    … …  
    284314      tlist.AddToList(&contdisp);
    285315    tlist.AddToList(&seltrain);
    286     tlist.AddToList(&filltrain);
     316    tlist.AddToList(&filltraincalc);
     317    tlist.AddToList(&filltrainhists);
    287318
    288319    //******************************
    … …  
    302333
    303334
    304     // print the filled Training Matrix
    305     fMatrixTrain->Print("SizeCols");
     335    // print the filled Training Matrices
     336    fMatrixTrainCalc->Print("SizeCols");
     337    fMatrixTrainHists->Print("SizeCols");
    306338
    307339    // check that number of generated events is compatible with the resquested
    308     Int_t howmanygenerated = fMatrixTrain->GetM().GetNrows();
    309     if (TMath::Abs(howmanygenerated-howmanytrain) > TMath::Sqrt(9.*howmanytrain))
     340    Int_t howmanygeneratedcalc = fMatrixTrainCalc->GetM().GetNrows();
     341    if (TMath::Abs(howmanygeneratedcalc-howmanytrain) > TMath::Sqrt(9.*howmanytrain))
    310342    {
    311343      *fLog << "MFindDisp::DefineTrainMatrix; no.of generated events ("
    312             << howmanygenerated
     344            << howmanygeneratedcalc
    313345            << ") is incompatible with the no.of requested events ("
    314346            << howmanytrain << ")" << endl;
    315347    }
    316348
    317     *fLog << "TRAINING Matrix was filled" << endl;
     349    Int_t howmanygeneratedhists = fMatrixTrainHists->GetM().GetNrows();
     350    if (TMath::Abs(howmanygeneratedhists-howmanytrain) > TMath::Sqrt(9.*howmanytrain))
     351    {
     352      *fLog << "MFindDisp::DefineTrainMatrix; no.of generated events ("
     353            << howmanygeneratedhists
     354            << ") is incompatible with the no.of requested events ("
     355            << howmanytrain << ")" << endl;
     356    }
     357
     358
     359    *fLog << "TRAINING Matrices were filled" << endl;
    318360    *fLog << "=============================================" << endl;
    319361
    320362
    321363    //--------------------------
    322     // write out training matrix
     364    // write out training matrices
    323365
    324366    if (filetrain != "")
    325367    {
    326368      TFile filetr(filetrain, "RECREATE");
    327       fMatrixTrain->Write();
     369      fMatrixTrainCalc->Write();
     370      fMatrixTrainHists->Write();
    328371      filetr.Close();
    329372
    330       *fLog << "MFindDisp::DefineTrainMatrix; Training matrix was written onto file '"
     373      *fLog << "MFindDisp::DefineTrainMatrix; Training matrices were written onto file '"
    331374            << filetrain << "'" << endl;
    332375    }
     376
    333377
    334378    if (display != NULL)
    … …  
    341385// --------------------------------------------------------------------------
    342386//
    343 // Define the matrix 'fMatrixTest' for the TEST sample
     387// Define the matrices 'fMatrixTestCalc' and 'fMatrixTestHists'
     388// for the TEST sample
    344389//
    345390// alltogether 'howmanytest' events are read from file 'nametest'
    … …  
    356401
    357402    *fLog   << "=============================================" << endl;
    358     *fLog   << "Fill TEST Matrix from file '" << nametest << endl;
     403    *fLog   << "Fill TEST Matrices from file '" << nametest << endl;
    359404    *fLog   << "     select " << howmanytest << " events " << endl;
    360405    if (!hreftest.GetHist().GetEntries()==0)
    … …  
    390435    seltest.SetName("selectTest");
    391436 
    392     MFillH filltest(fMatrixTest);
    393     filltest.SetFilter(&seltest);
    394     filltest.SetName("fillMatrixTest");
     437    MFillH filltestcalc(fMatrixTestCalc);
     438    filltestcalc.SetFilter(&seltest);
     439    filltestcalc.SetName("fillMatrixTestCalc");
     440
     441    MFillH filltesthists(fMatrixTestHists);
     442    filltesthists.SetFilter(&seltest);
     443    filltesthists.SetName("fillMatrixTestHists");
    395444
    396445    //******************************
    … …  
    399448    plist.AddToList(&tlist);
    400449    plist.AddToList(fCam);
    401     plist.AddToList(fMatrixTest);
     450    plist.AddToList(fMatrixTestCalc);
     451    plist.AddToList(fMatrixTestHists);
    402452
    403453    //******************************
    … …  
    408458      tlist.AddToList(&contdisp);
    409459    tlist.AddToList(&seltest);
    410     tlist.AddToList(&filltest);
     460    tlist.AddToList(&filltestcalc);
     461    tlist.AddToList(&filltesthists);
    411462
    412463    //******************************
    … …  
    426477
    427478
    428     // print the filled Test Matrix
    429     fMatrixTest->Print("SizeCols");
     479    // print the filled Test Matrices
     480    fMatrixTestCalc->Print("SizeCols");
     481    fMatrixTestHists->Print("SizeCols");
    430482
    431483    // check that number of generated events is compatible with the resquested
    432     const Int_t howmanygenerated = fMatrixTest->GetM().GetNrows();
    433     if (TMath::Abs(howmanygenerated-howmanytest) > TMath::Sqrt(9.*howmanytest))
     484    const Int_t howmanygeneratedcalc = fMatrixTestCalc->GetM().GetNrows();
     485    if (TMath::Abs(howmanygeneratedcalc-howmanytest) > TMath::Sqrt(9.*howmanytest))
    434486    {
    435487      *fLog << "MFindDisp::DefineTestMatrix; no.of generated events ("
    436             << howmanygenerated 
     488            << howmanygeneratedcalc
    437489            << ") is incompatible with the no.of requested events ("
    438490            << howmanytest << ")" << endl;
    439491    }
    440492
    441     *fLog << "TEST Matrix was filled" << endl;
     493    const Int_t howmanygeneratedhists = fMatrixTestHists->GetM().GetNrows();
     494    if (TMath::Abs(howmanygeneratedhists-howmanytest) > TMath::Sqrt(9.*howmanytest))
     495    {
     496      *fLog << "MFindDisp::DefineTestMatrix; no.of generated events ("
     497            << howmanygeneratedhists
     498            << ") is incompatible with the no.of requested events ("
     499            << howmanytest << ")" << endl;
     500    }
     501
     502    *fLog << "TEST Matrices were filled" << endl;
    442503    *fLog << "=============================================" << endl;
    443504
    444505
    445506    //--------------------------
    446     // write out test matrix
     507    // write out test matrices
    447508
    448509    if (filetest != "")
    449510    {
    450511      TFile filete(filetest, "RECREATE");
    451       fMatrixTest->Write();
     512      fMatrixTestCalc->Write();
     513      fMatrixTestHists->Write();
    452514      filete.Close();
    453515
    454       *fLog << "MFindDisp::DefineTestMatrix; Test matrix was written onto file '"
     516      *fLog << "MFindDisp::DefineTestMatrix; Test matrices were written onto file '"
    455517            << filetest << "'" << endl;
    456518    }
     519
    457520
    458521    if (display != NULL)
    … …  
    477540{
    478541    *fLog   << "=============================================" << endl;
    479     *fLog   << "Fill TRAINING and TEST Matrix from file '" << name << endl;
     542    *fLog   << "Fill TRAINING and TEST Matrices from file '" << name << endl;
    480543    *fLog   << "     select "   << howmanytrain
    481544            << " training and " << howmanytest << " test events " << endl;
    … …  
    521584    split.SetSelectionRatio(prob);
    522585
    523     MFillH filltrain(fMatrixTrain);
    524     filltrain.SetFilter(&split);
    525     filltrain.SetName("fillMatrixTrain");
     586    MFillH filltraincalc(fMatrixTrainCalc);
     587    filltraincalc.SetFilter(&split);
     588    filltraincalc.SetName("fillMatrixTrainCalc");
     589
     590    MFillH filltrainhists(fMatrixTrainHists);
     591    filltrainhists.SetFilter(&split);
     592    filltrainhists.SetName("fillMatrixTrainHists");
    526593
    527594    // consider this event as candidate for a test event
    … …  
    530597    conttrain.SetName("ContTrain");
    531598
    532     MFillH filltest(fMatrixTest);
    533     filltest.SetName("fillMatrixTest");
     599    MFillH filltestcalc(fMatrixTestCalc);
     600    filltestcalc.SetName("fillMatrixTestCalc");
     601
     602    MFillH filltesthists(fMatrixTestHists);
     603    filltesthists.SetName("fillMatrixTestHists");
    534604
    535605
    … …  
    539609    plist.AddToList(&tlist);
    540610    plist.AddToList(fCam);
    541     plist.AddToList(fMatrixTrain);
    542     plist.AddToList(fMatrixTest);
     611    plist.AddToList(fMatrixTrainCalc);
     612    plist.AddToList(fMatrixTrainHists);
     613    plist.AddToList(fMatrixTestCalc);
     614    plist.AddToList(fMatrixTestHists);
    543615
    544616    //******************************
    … …  
    549621      tlist.AddToList(&contdisp);
    550622    tlist.AddToList(&cont);
    551     tlist.AddToList(&filltrain);
     623    tlist.AddToList(&filltraincalc);
     624    tlist.AddToList(&filltrainhists);
    552625    tlist.AddToList(&conttrain);
    553     tlist.AddToList(&filltest);
     626    tlist.AddToList(&filltestcalc);
     627    tlist.AddToList(&filltesthists);
    554628
    555629    //******************************
    … …  
    570644
    571645
    572     // print the filled Train Matrix
    573     fMatrixTrain->Print("SizeCols");
     646    // print the filled Train Matrices
     647    fMatrixTrainCalc->Print("SizeCols");
     648    fMatrixTrainHists->Print("SizeCols");
    574649
    575650    // check that number of generated events is compatible with the resquested
    576     const Int_t generatedtrain = fMatrixTrain->GetM().GetNrows();
    577     if (TMath::Abs(generatedtrain-howmanytrain) > TMath::Sqrt(9.*howmanytrain))
     651    const Int_t generatedtraincalc = fMatrixTrainCalc->GetM().GetNrows();
     652    if (TMath::Abs(generatedtraincalc-howmanytrain) > TMath::Sqrt(9.*howmanytrain))
    578653    {
    579654      *fLog << "MFindDisp::DefineTrainTestMatrix; no.of generated events ("
    580             << generatedtrain 
     655            << generatedtraincalc
    581656            << ") is incompatible with the no.of requested events ("
    582657            << howmanytrain << ")" << endl;
    583658    }
    584 
    585 
    586     // print the filled Test Matrix
    587     fMatrixTest->Print("SizeCols");
     659    const Int_t generatedtrainhists = fMatrixTrainHists->GetM().GetNrows();
     660    if (TMath::Abs(generatedtrainhists-howmanytrain) > TMath::Sqrt(9.*howmanytrain))
     661    {
     662      *fLog << "MFindDisp::DefineTrainTestMatrix; no.of generated events ("
     663            << generatedtrainhists
     664            << ") is incompatible with the no.of requested events ("
     665            << howmanytrain << ")" << endl;
     666    }
     667
     668
     669    // print the filled Test Matrices
     670    fMatrixTestCalc->Print("SizeCols");
     671    fMatrixTestHists->Print("SizeCols");
    588672
    589673    // check that number of generated events is compatible with the resquested
    590     const Int_t generatedtest = fMatrixTest->GetM().GetNrows();
    591     if (TMath::Abs(generatedtest-howmanytest) > TMath::Sqrt(9.*howmanytest))
     674    const Int_t generatedtestcalc = fMatrixTestCalc->GetM().GetNrows();
     675    if (TMath::Abs(generatedtestcalc-howmanytest) > TMath::Sqrt(9.*howmanytest))
    592676    {
    593677      *fLog << "MFindDisp::DefineTrainTestMatrix; no.of generated events ("
    594             << generatedtest 
     678            << generatedtestcalc
    595679            << ") is incompatible with the no.of requested events ("
    596680            << howmanytest << ")" << endl;
    597681    }
    598 
    599 
    600     *fLog << "TRAINING and TEST Matrix were filled" << endl;
     682    const Int_t generatedtesthists = fMatrixTestHists->GetM().GetNrows();
     683    if (TMath::Abs(generatedtesthists-howmanytest) > TMath::Sqrt(9.*howmanytest))
     684    {
     685      *fLog << "MFindDisp::DefineTrainTestMatrix; no.of generated events ("
     686            << generatedtesthists
     687            << ") is incompatible with the no.of requested events ("
     688            << howmanytest << ")" << endl;
     689    }
     690
     691
     692    *fLog << "TRAINING and TEST Matrices were filled" << endl;
    601693    *fLog << "=============================================" << endl;
    602694
    603695
     696    //----------------------------
     697    // write out training matrices
     698
     699    if (filetrain != "")
     700    {
     701      TFile filetr(filetrain, "RECREATE");
     702      fMatrixTrainCalc->Write();
     703      fMatrixTrainHists->Write();
     704      filetr.Close();
     705
     706      *fLog << "MFindDisp::DefineTrainTestMatrix; Training matrices were written onto file '"
     707            << filetrain << "'" << endl;
     708    }
     709
    604710    //--------------------------
    605     // write out training matrix
    606 
    607     if (filetrain != "")
    608     {
    609       TFile filetr(filetrain, "RECREATE");
    610       fMatrixTrain->Write();
    611       filetr.Close();
    612 
    613       *fLog << "MFindDisp::DefineTrainTestMatrix; Training matrix was written onto file '"
    614             << filetrain << "'" << endl;
    615     }
    616 
    617     //--------------------------
    618     // write out test matrix
     711    // write out test matrices
    619712
    620713    if (filetest != "")
    621714    {
    622715      TFile filete(filetest, "RECREATE");
    623       fMatrixTest->Write();
     716      fMatrixTestCalc->Write();
     717      fMatrixTestHists->Write();
    624718      filete.Close();
    625719
    626       *fLog << "MFindDisp::DefineTrainTestMatrix; Test matrix was written onto file '"
     720      *fLog << "MFindDisp::DefineTrainTestMatrix; Test matrices were written onto file '"
    627721            << filetest << "'" << endl;
    628722    }
    … …  
    643737{
    644738  //--------------------------
    645   // read in training matrix
     739  // read in training matrices
    646740
    647741  TFile filetr(filetrain);
    648   fMatrixTrain->Read("MatrixTrain");
    649   fMatrixTrain->Print("SizeCols");
    650 
    651   *fLog << "MFindDisp::ReadMatrix; Training matrix was read in from file '"
     742  fMatrixTrainCalc->Read("MatrixTrainCalc");
     743  fMatrixTrainHists->Read("MatrixTrainHists");
     744  fMatrixTrainCalc->Print("SizeCols");
     745  fMatrixTrainHists->Print("SizeCols");
     746
     747  *fLog << "MFindDisp::ReadMatrix; Training matrices were read in from file '"
    652748        << filetrain << "'" << endl;
    653749  filetr.Close();
    … …  
    655751
    656752  //--------------------------
    657   // read in test matrix
     753  // read in test matrices
    658754
    659755  TFile filete(filetest);
    660   fMatrixTest->Read("MatrixTest");
    661   fMatrixTest->Print("SizeCols");
    662 
    663   *fLog << "MFindDisp::ReadMatrix; Test matrix was read in from file '"
     756  fMatrixTestCalc->Read("MatrixTestCalc");
     757  fMatrixTestHists->Read("MatrixTestHists");
     758  fMatrixTestCalc->Print("SizeCols");
     759  fMatrixTestHists->Print("SizeCols");
     760
     761  *fLog << "MFindDisp::ReadMatrix; Test matrices were read in from file '"
    664762        << filetest << "'" << endl;
    665763  filete.Close();
    … …  
    686784//
    687785// - call TMinuit to do the minimization :
    688 //        the fcnDisp function calculates the Chi^2 for the current
    689 //                                                parameter values;
     786//        the fcnDisp function calculates the parameter to minimize
     787//                            for the current Disp parameter values;
    690788//        for this - Disp is calculated in the event loop by calling
    691 //                   MDispCalc::Process() ==> MDisp::Calc()
    692 //                 - the Chi^2 contributions are summed up in the event loop
    693 //                   by calling MHDisp::Fill()
    694 //                 - after the event loop the Chi^2 is calculated by calling
    695 //                   MHDisp::Finalize()
     789//                   MDispCalc::Process() ==> MDispCalc::Calc()
     790//                 - the Minimization parameter contributions are summed up
     791//                   in the event loop by calling MHDisp::Fill()
     792//                 - after the event loop the final value of the Minimization
     793//                   parameter is calculated by calling MHDisp::Finalize()
    696794//
    697795// Needed as input : (to be set by the Set functions)
    … …  
    703801//     - from the file parDispInit (if it is != "")
    704802//     - or from the arrays params and/or steps
    705 //     - or from the Disp constructor
     803//     - or from the DispParameters constructor
    706804//
    707805//----------------------------------------------------------------------
    … …  
    721819
    722820
    723     if (fMatrixTrain == NULL)
    724     {
    725       *fLog << "MFindDisp::FindParams; training matrix is not defined... aborting"
     821    if (fMatrixTrainCalc == NULL || fMatrixTrainHists == NULL)
     822    {
     823      *fLog << "MFindDisp::FindParams; training matrices are not defined... aborting"
    726824            << endl;
    727825      return kFALSE;
    728826    }
    729827
    730     if (fMatrixTrain->GetM().GetNrows() <= 0)
    731     {
    732       *fLog << "MFindDisp::FindParams; training matrix has no entries"
     828    if (fMatrixTrainCalc->GetM().GetNrows() <= 0  ||  fMatrixTrainHists->GetM().GetNrows() <= 0)
     829    {
     830      *fLog << "MFindDisp::FindParams; training matrices have no entries"
    733831            << endl;
    734832      return kFALSE;
    … …  
    740838
    741839    // loop over rows of matrix
    742     MMatrixLoop loop(fMatrixTrain);
     840    MMatrixLoop loopcalc(fMatrixTrainCalc);
     841    MMatrixLoop loophists(fMatrixTrainHists);
    743842
    744843    //--------------------------------
    745844    // create container for the Disp parameters
    746845    // and set them to their initial values
    747     MDisp disp;
     846    MDispParameters dispparams;
    748847
    749848    // take initial values from file parDispInit
    … …  
    751850    {
    752851      TFile inparam(parDispInit);
    753       disp.Read("MDisp");
     852      dispparams.Read("MDispParameters");
    754853      inparam.Close();
    755854      *fLog << "MFindDisp::FindParams; initial values of parameters are taken from file "
    … …  
    764863        *fLog << "MFindDisp::FindParams; initial values of parameters are taken from 'params'"
    765864              << endl;
    766         disp.SetParameters(params);
     865        dispparams.SetParameters(params);
    767866      }
    768867      if (steps.GetSize()  != 0)
    … …  
    770869        *fLog << "MFindDisp::FindParams; initial step sizes are taken from 'steps'"
    771870              << endl;
    772         disp.SetStepsizes(steps);
     871        dispparams.SetStepsizes(steps);
    773872      }
    774873    }
    … …  
    776875    {
    777876        *fLog << "MFindDisp::FindParams; initial values and step sizes are taken "
    778               << "from the MDisp constructor" << endl;
    779     }
    780 
    781     // get the matrix pointer and mapping from MDispCalc object to set the same in MHDisp object
    782     TArrayI map;
    783     MHMatrix *tmpmatrix = NULL;
    784     tmpmatrix = fDispCalcTrain->GetMatrixMap(map);
    785 
    786     // attention: argument of MHDisp is the name of MImageParDisp container, that should
    787     // be the same than the name given to it when creating MDispCalc object at the MFindDisp
    788     // constructor:  fDispCalcTrain = new MDispCalc("DispTrain");
    789     MHDisp hdisp("DispTrain");
    790     hdisp.SetMatrixMap(tmpmatrix,map);
    791     // fill the plots for Disp and sum up the Chi^2 contributions
     877              << "from the MDispParameters constructor" << endl;
     878    }
     879
     880    // fill the plots for Disp and sum up the Minimization parameter contributions
    792881    MFillH filldispplots("MHDisp", "");
    793882
    … …  
    796885   
    797886    parlistfcn.AddToList(&tasklistfcn);
    798     parlistfcn.AddToList(&disp);
    799     parlistfcn.AddToList(&hdisp);
     887    parlistfcn.AddToList(&dispparams);
     888    parlistfcn.AddToList(fHDispTrain);
    800889    parlistfcn.AddToList(fCam);
    801     parlistfcn.AddToList(fMatrixTrain);
     890    parlistfcn.AddToList(fMatrixTrainCalc);
     891    parlistfcn.AddToList(fMatrixTrainHists);
    802892
    803893    //******************************
    804894    // entries in MTaskList
    805895
    806     tasklistfcn.AddToList(&loop);
     896    tasklistfcn.AddToList(&loopcalc);
     897    tasklistfcn.AddToList(&loophists);
    807898    tasklistfcn.AddToList(fDispCalcTrain);
    808899    tasklistfcn.AddToList(&filldispplots);
    … …  
    834925
    835926    // get initial values of parameters
    836     fVinit = disp.GetParameters();
    837     fStep  = disp.GetStepsizes();
     927    fVinit = dispparams.GetParameters();
     928    fStep  = dispparams.GetStepsizes();
    838929
    839930    TString name[fVinit.GetSize()];
    … …  
    901992   
    902993    TFile outparam(fFilenameParam, "RECREATE");
    903     disp.Write();
     994    dispparams.Write();
    904995    outparam.Close();
    905996
    … …  
    907998              << fFilenameParam << "' :" << endl;
    908999
    909     const TArrayD &check = disp.GetParameters();
     1000    const TArrayD &check = dispparams.GetParameters();
    9101001    for (Int_t i=0; i<check.GetSize(); i++)
    9111002        *fLog << check[i] << ",  ";
    … …  
    9161007    //-------------------------------------------
    9171008    // draw the plots
    918     hdisp.Draw();
     1009    fHDispTrain->Draw();
    9191010
    9201011    *fLog << "End of Disp optimization part" << endl;
    … …  
    9321023Bool_t MFindDisp::TestParams()
    9331024{
    934     if (fMatrixTest == NULL)
    935     {
    936       *fLog << "MFindDisp::TestParams; test matrix is not defined... aborting"
     1025    if (fMatrixTestCalc == NULL || fMatrixTestHists == NULL)
     1026    {
     1027      *fLog << "MFindDisp::TestParams; test matrices are not defined... aborting"
    9371028            << endl;
    9381029      return kFALSE;
    9391030    }
    9401031
    941     if (fMatrixTest->GetM().GetNrows() <= 0)
    942     {
    943         *fLog << "MFindDisp::TestParams; test matrix has no entries"
     1032    if (fMatrixTestCalc->GetM().GetNrows() <= 0  ||  fMatrixTestHists->GetM().GetNrows() <= 0)
     1033    {
     1034        *fLog << "MFindDisp::TestParams; test matrices have no entries"
    9441035              << endl;
    9451036        return kFALSE;
    … …  
    9551046
    9561047    TFile inparam(fFilenameParam);
    957     MDisp dispin;
    958     dispin.Read("MDisp");
     1048    MDispParameters dispin;
     1049    dispin.Read("MDispParameters");
    9591050    inparam.Close();
    9601051
    … …  
    9731064    MTaskList tasklist2;
    9741065
    975     MDisp disp;
    976     disp.SetParameters(dispPar);
    977 
    978     MMatrixLoop loop(fMatrixTest);
    979 
    980     // get the matrix pointer and mapping from MDispCalc object to set the same in MHDisp object
    981     TArrayI map;
    982     MHMatrix *tmpmatrix = NULL;
    983     tmpmatrix = fDispCalcTest->GetMatrixMap(map);
     1066    MDispParameters dispparams;
     1067    dispparams.SetParameters(dispPar);
     1068
     1069    MMatrixLoop loopcalc(fMatrixTestCalc);
     1070    MMatrixLoop loophists(fMatrixTestHists);
     1071
     1072    MHMatrix *imatrix = NULL;
     1073    TArrayI  imap;
     1074    fHDispTest->GetMatrixMap(imatrix,imap);
    9841075
    9851076    // attention: argument of MHDisp is the name of MImageParDisp container, that should
    9861077    // be the same than the name given to it when creating MDispCalc object at the MFindDisp
    9871078    // constructor:  fDispCalcTrain = new MDispCalc("DispTest");
    988     // fill the plots for Disp and sum up the Chi^2 contributions
     1079    // fill the plots for Disp and sum up the Minimization parameter contributions
    9891080    MHDisp hdisp1("DispTest");
    9901081    hdisp1.SetName("MHDispCorr");
    9911082    hdisp1.SetSelectedPos(1);
    992     hdisp1.SetMatrixMap(tmpmatrix,map);
     1083    hdisp1.SetMatrixMap(imatrix,imap);
    9931084    MFillH filldispplots1("MHDispCorr[MHDisp]", "");
    9941085   
    … …  
    9961087    hdisp2.SetName("MHDispWrong");
    9971088    hdisp2.SetSelectedPos(2);
    998     hdisp2.SetMatrixMap(tmpmatrix,map);
     1089    hdisp2.SetMatrixMap(imatrix,imap);
    9991090    MFillH filldispplots2("MHDispWrong[MHDisp]", "");
    10001091   
    … …  
    10021093    hdisp3.SetName("MHDispM3Long");
    10031094    hdisp3.SetSelectedPos(3);
    1004     hdisp3.SetMatrixMap(tmpmatrix,map);
     1095    hdisp3.SetMatrixMap(imatrix,imap);
    10051096    MFillH filldispplots3("MHDispM3Long[MHDisp]", "");
    10061097   
    … …  
    10081099    hdisp4.SetName("MHDispAsym");
    10091100    hdisp4.SetSelectedPos(4);
    1010     hdisp4.SetMatrixMap(tmpmatrix,map);
     1101    hdisp4.SetMatrixMap(imatrix,imap);
    10111102    MFillH filldispplots4("MHDispAsym[MHDisp]", "");
    10121103   
    … …  
    10161107
    10171108    parlist2.AddToList(&tasklist2);
    1018     parlist2.AddToList(&disp);
     1109    parlist2.AddToList(&dispparams);
    10191110    parlist2.AddToList(&hdisp1);
    10201111    parlist2.AddToList(&hdisp2);
    … …  
    10221113    parlist2.AddToList(&hdisp4);
    10231114    parlist2.AddToList(fCam);
    1024     parlist2.AddToList(fMatrixTest);
     1115    parlist2.AddToList(fMatrixTestCalc);
     1116    parlist2.AddToList(fMatrixTestHists);
    10251117
    10261118    //******************************
    10271119    // entries in MTaskList
    10281120
    1029     tasklist2.AddToList(&loop);
     1121    tasklist2.AddToList(&loopcalc);
     1122    tasklist2.AddToList(&loophists);
    10301123    tasklist2.AddToList(fDispCalcTest);
    10311124    tasklist2.AddToList(&filldispplots1);
    … …  
    10641157    return kTRUE;
    10651158}
    1066 
    1067 
    1068 
    1069 
    1070 
    1071 
    1072 
    1073 

  • trunk/MagicSoft/Mars/mtemp/mifae/library/MFindDisp.h

    r5671 r5904  
    1616class MHMatrix;
    1717class MDispCalc;
     18class MHDisp;
    1819class MFDisp;
    1920class MFilter;
    … …  
    3435  Int_t     fHowManyTest;    // number of events for testing
    3536
    36   MHMatrix  *fMatrixTrain;   // training matrix
    37   MHMatrix  *fMatrixTest;    // testing matrix
     37  MHMatrix  *fMatrixTrainCalc;    // training matrix with variables needed in MDispCalc
     38  MHMatrix  *fMatrixTrainHists;   // training matrix with variables needed in MHDisp
     39  MHMatrix  *fMatrixTestCalc;     // testing matrix with variables needed in MDispCalc
     40  MHMatrix  *fMatrixTestHists;    // testing matrix with variables needed in MHDisp
    3841
    3942  MDispCalc *fDispCalcTrain;
    4043  MDispCalc *fDispCalcTest;
     44
     45  MHDisp    *fHDispTrain;
     46  MHDisp    *fHDispTest;
    4147
    4248  MFDisp    *fDispFilter;    // filter to select an events sample
    … …  
    9399                               Int_t iflag=0);
    94100
    95   MHMatrix *GetMatrixTrain() { return fMatrixTrain; }
    96   MHMatrix *GetMatrixTest()  { return fMatrixTest;  }
     101  MHMatrix *GetMatrixTrainCalc()  { return fMatrixTrainCalc;  }
     102  MHMatrix *GetMatrixTrainHists() { return fMatrixTrainHists; }
     103  MHMatrix *GetMatrixTestCalc()   { return fMatrixTestCalc;   }
     104  MHMatrix *GetMatrixTestHists()  { return fMatrixTestHists;  }
    97105
    98106  Bool_t ReadMatrix( const TString &filetrain, const TString &filetest);

  • trunk/MagicSoft/Mars/mtemp/mifae/library/MHDisp.cc

    r5671 r5904  
    2929//                                                                         //
    3030//   container holding the histograms for Disp                             //
    31 //   also computes the Chi^2 of the Disp optimization                      //
     31//   also computes the minimization parameter of the Disp optimization     //
    3232//                                                                         //
    3333/////////////////////////////////////////////////////////////////////////////
    … …  
    5454
    5555#include "MHMatrix.h"
     56#include "MParameters.h"
    5657
    5758#include "MParList.h"
    … …  
    6364
    6465using namespace std;
     66
     67enum dispVar_t {kX,kY,kMeanX,kMeanY,kDelta,kSize,kM3Long,kAsym,
     68                kEnergy,kImpact,kLongitmax,kZd,kAz,kTotVar};    // enum variables for the
     69                                                                // matrix columns mapping
    6570
    6671// --------------------------------------------------------------------------
    … …  
    7580    fTitle = title ? title : "Histograms for Disp";
    7681
    77     fSelectedPos = 1;  // default MHDisp flag (selects Correct Disp srcpos)
     82    fSelectedPos = 1;  // default MHDisp flag (selects Correct Disp source position solution)
    7883
    7984    fMatrix = NULL;
     85
     86    // initialize mapping array dimension to the number of columns of fMatrix
     87    fMap.Set(kTotVar);
    8088
    8189    //--------------------------------------------------
    … …  
    110118    fSkymapXY->SetYTitle("Y Disp source position [deg]");
    111119
    112     fHistChi2   = new TH1F("fHistChi2"  ,
     120    fHistMinPar   = new TH1F("fHistMinPar"  ,
    113121         "Distance^2 between Disp and real srcpos", 100, 0., 2.);
    114     fHistChi2->SetDirectory(NULL);
    115     fHistChi2->UseCurrentStyle();
    116     fHistChi2->SetXTitle("Chi^2 = distance^2 Disp to real srcpos [deg^2]");
    117     fHistChi2->SetYTitle("# events");
     122    fHistMinPar->SetDirectory(NULL);
     123    fHistMinPar->UseCurrentStyle();
     124    fHistMinPar->SetXTitle("Minimization parameter = d^2 Disp to real srcpos [deg^2]");
     125    fHistMinPar->SetYTitle("# events");
    118126   
    119127    fHistDuDv   = new TH2F("fHistDuDv",
    … …  
    125133    fHistDuDv->SetYTitle("Du = longitudinal distance [deg]");
    126134
    127     fHistChi2Energy   = new TH2F("fHistChi2Energy",
    128          "Chi^2 vs. Energy", 50, 1., 3., 100, 0., 2.);
    129     fHistChi2Energy->SetDirectory(NULL);
    130     fHistChi2Energy->UseCurrentStyle();
    131     fHistChi2Energy->SetXTitle("log10(Energy) [GeV]");
    132     fHistChi2Energy->SetYTitle("Chi^2 = distance^2 Disp to real srcpos [deg^2]");
    133 
    134     fHistChi2Size   = new TH2F("fHistChi2Size",
    135          "Chi^2 vs. Size", 50, 2., 4., 100, 0., 2.);
    136     fHistChi2Size->SetDirectory(NULL);
    137     fHistChi2Size->UseCurrentStyle();
    138     fHistChi2Size->SetXTitle("log10(Size) [#phot]");
    139     fHistChi2Size->SetYTitle("Chi^2 = distance^2 Disp to real srcpos [deg^2]");
     135    fHistMinParEnergy   = new TH2F("fHistMinParEnergy",
     136         "Minimization parameter vs. Energy", 50, 1., 3., 100, 0., 2.);
     137    fHistMinParEnergy->SetDirectory(NULL);
     138    fHistMinParEnergy->UseCurrentStyle();
     139    fHistMinParEnergy->SetXTitle("log10(Energy) [GeV]");
     140    fHistMinParEnergy->SetYTitle("Minimization parameter = d^2 Disp to real srcpos [deg^2]");
     141
     142    fHistMinParSize   = new TH2F("fHistMinParSize",
     143         "Minimization parameter vs. Size", 50, 2., 4., 100, 0., 2.);
     144    fHistMinParSize->SetDirectory(NULL);
     145    fHistMinParSize->UseCurrentStyle();
     146    fHistMinParSize->SetXTitle("log10(Size) [#phot]");
     147    fHistMinParSize->SetYTitle("Minimization parameter = d^2 Disp to real srcpos [deg^2]");
    140148
    141149    fHistDuEnergy   = new TH2F("fHistDuEnergy",
    … …  
    170178         "Fraction events srcpos well assign vs. Size", 50, 2., 4., 0., 1.);
    171179    fEvCorrAssign->SetDirectory(NULL);
    172     fEvCorrAssign->UseCurrentStyle();
     180    fEvCorrAssign->SetStats(0);
    173181    fEvCorrAssign->SetXTitle("log10(Size) [#phot]");
    174182    fEvCorrAssign->SetYTitle("Fraction of events with srcpos well assign");
    … …  
    186194  delete fHistcosZA;
    187195  delete fSkymapXY;
    188   delete fHistChi2;
     196  delete fHistMinPar;
    189197  delete fHistDuDv;
    190   delete fHistChi2Energy;
    191   delete fHistChi2Size;
     198  delete fHistMinParEnergy;
     199  delete fHistMinParSize;
    192200  delete fHistDuEnergy;
    193201  delete fHistDuSize;
    … …  
    204212Bool_t MHDisp::SetupFill(const MParList *pList)
    205213{
    206     // reset all histograms and Chi^2 computing variables
     214    // reset all histograms and Minimization parameter computing variables
    207215    // before each new eventloop
    208216    fNumEv = 0;
    209     fSumChi2  = 0.;
    210     fChi2  = 0.;
     217    fSumMinPar  = 0.;
     218    fMinPar = (MParameterD*)const_cast<MParList*>(pList)->FindCreateObj("MParameterD", "MinimizationParameter");
     219    if (!fMinPar)
     220    {
     221      *fLog << err << "MParameterD (MinimizationParameter) not found and could not be created... aborting."
     222            << endl;
     223      return kFALSE;
     224    }
     225    fMinPar->SetVal(0);
    211226
    212227    fHistEnergy->Reset();
    … …  
    214229    fHistcosZA->Reset();
    215230    fSkymapXY->Reset();
    216     fHistChi2->Reset();
     231    fHistMinPar->Reset();
    217232    fHistDuDv->Reset();
    218     fHistChi2Energy->Reset();
    219     fHistChi2Size->Reset();
     233    fHistMinParEnergy->Reset();
     234    fHistMinParSize->Reset();
    220235    fHistDuEnergy->Reset();
    221236    fHistDuSize->Reset();
    … …  
    250265    // if fMatrix exists, skip preprocessing and just read events from matrix;
    251266    // if doesn't exist, read variables values from containers, so look for them
     267    cout << "MHDisp::SetupFill:  fMatrix = " << fMatrix << endl;
     268
    252269    if (fMatrix)
    253270      return kTRUE;
    … …  
    291308    if (!fPointing)
    292309    {
    293         *fLog << err << "MPointingPos not found... "
    294               << "MMcEvt is going to be used to get Theta and Phi."
    295               << endl;
    296         //        return kFALSE;
     310        *fLog << err << "MPointingPos not found... aborting." << endl;
     311        return kFALSE;
    297312    }
    298313
    … …  
    306321//
    307322// Set which selection algorithm for the Disp estimated source position
    308 // we want to follow when filling the histograms:
     323// solutions we want to follow when filling the histograms:
    309324//  * iflag == 1 --> Correct source position, the closest to the real srcpos
    310325//                   (only applicable to MC or well known source real data)
    … …  
    323338//
    324339// Sets the Matrix and the array of mapped values for each Matrix column
    325 // (see MDispCalc)
    326340//
    327341void MHDisp::SetMatrixMap(MHMatrix *matrix, TArrayI &map)
    … …  
    334348// --------------------------------------------------------------------------
    335349//
     350// Returns the Matrix and the mapped value for each Matrix column
     351//
     352void MHDisp::GetMatrixMap(MHMatrix* &matrix, TArrayI &map)
     353{
     354    map = fMap;
     355    matrix = fMatrix;
     356}
     357
     358
     359// --------------------------------------------------------------------------
     360//
     361// You can use this function if you want to use a MHMatrix instead of the
     362// given containers for the Disp optimization. This function adds all
     363// necessary columns to the given matrix. Afterwards, you should fill
     364// the matrix with the corresponding data (eg from a file by using
     365// MHMatrix::Fill). Then, if you loop through the matrix (eg using
     366// MMatrixLoop), MHDisp::Fill will take the values from the matrix
     367// instead of the containers.
     368//
     369void MHDisp::InitMapping(MHMatrix *mat)
     370{
     371    if (fMatrix)
     372      return;
     373
     374    fMatrix = mat;
     375
     376    fMap[kX]          = fMatrix->AddColumn("MSrcPosCam.fX");
     377    fMap[kY]          = fMatrix->AddColumn("MSrcPosCam.fY");
     378    fMap[kMeanX]      = fMatrix->AddColumn("MHillas.fMeanX");
     379    fMap[kMeanY]      = fMatrix->AddColumn("MHillas.fMeanY");
     380    fMap[kDelta]      = fMatrix->AddColumn("MHillas.fDelta");
     381
     382    fMap[kSize]       = fMatrix->AddColumn("MHillas.fSize");
     383   
     384    fMap[kM3Long]     = fMatrix->AddColumn("MHillasExt.fM3Long");
     385    fMap[kAsym]       = fMatrix->AddColumn("MHillasExt.fAsym");
     386   
     387    fMap[kEnergy]     = fMatrix->AddColumn("MMcEvt.fEnergy");
     388    fMap[kImpact]     = fMatrix->AddColumn("MMcEvt.fImpact");
     389    fMap[kLongitmax]  = fMatrix->AddColumn("MMcEvt.fLongitmax");
     390   
     391    fMap[kZd]         = fMatrix->AddColumn("MPointingPos.fZd");
     392    fMap[kAz]         = fMatrix->AddColumn("MPointingPos.fAz");
     393}
     394
     395
     396// --------------------------------------------------------------------------
     397//
    336398// Returns a mapped value from the Matrix
    337399//
    … …  
    349411// --------------------------------------------------------------------------
    350412//
    351 // Fill the histograms and sum up the Chi^2 outcome of each processed event
     413// Fill the histograms and sum up the Minimization paramter outcome
     414// of each processed event
    352415//
    353416Bool_t MHDisp::Fill(const MParContainer *par, const Stat_t w)
    … …  
    356419    Double_t impact = 0.;
    357420    Double_t xmax = 0.;
    358     Double_t theta = 0.;
    359     Double_t phi = 0.;
    360421
    361422    if ( fMatrix || (!fMatrix && fMcEvt) )
    362423    { 
    363       energy   = fMatrix ? GetVal(13)  : fMcEvt->GetEnergy();
    364       impact   = fMatrix ? GetVal(14)  : fMcEvt->GetImpact();
    365       xmax     = fMatrix ? GetVal(15)  : fMcEvt->GetLongitmax();
    366      
    367       if (fPointing)
    368       {
    369         theta  = fMatrix ? GetVal(16)  : fPointing->GetZd();
    370         phi    = fMatrix ? GetVal(17)  : fPointing->GetAz();
    371       }
    372       else
    373       {
    374         theta  = fMatrix ? GetVal(16)  : fMcEvt->GetTelescopeTheta()*kRad2Deg;
    375         phi    = fMatrix ? GetVal(17)  : fMcEvt->GetTelescopePhi()*kRad2Deg;
    376       }
     424      energy   = fMatrix ? GetVal(kEnergy)     : fMcEvt->GetEnergy();
     425      impact   = fMatrix ? GetVal(kImpact)     : fMcEvt->GetImpact();
     426      xmax     = fMatrix ? GetVal(kLongitmax)  : fMcEvt->GetLongitmax();     
    377427    }
    378428    else
    379429      *fLog << "No MMcEvt container available (no Energy,ImpactPar or Xmax)" << endl;
    380430
    381     Double_t xsrcpos0   = fMatrix ? GetVal(0)   : fSrcPos->GetX();
    382     Double_t ysrcpos0   = fMatrix ? GetVal(1)   : fSrcPos->GetY();
    383     Double_t xmean0     = fMatrix ? GetVal(2)   : fHil->GetMeanX(); 
    384     Double_t ymean0     = fMatrix ? GetVal(3)   : fHil->GetMeanY();
    385     Double_t delta      = fMatrix ? GetVal(4)   : fHil->GetDelta();
    386    
    387     Double_t size       = fMatrix ? GetVal(5)   : fHil->GetSize();
    388     //    Double_t width0     = fMatrix ? GetVal(6)   : fHil->GetWidth();
    389     //    Double_t length0    = fMatrix ? GetVal(7)   : fHil->GetLength();
    390    
    391     //    Double_t conc       = fMatrix ? GetVal(8)   : fNewPar->GetConc();
    392     //    Double_t leakage1   = fMatrix ? GetVal(9)   : fNewPar->GetLeakage1();
    393     //    Double_t leakage2   = fMatrix ? GetVal(10)  : fNewPar->GetLeakage2();
    394    
    395     Double_t m3long     = fMatrix ? GetVal(11)  : fHilExt->GetM3Long();
    396     Double_t asym       = fMatrix ? GetVal(12)  : fHilExt->GetAsym();
     431    Double_t theta      = fMatrix ? GetVal(kZd)  : fPointing->GetZd();
     432    //    Double_t phi        = fMatrix ? GetVal(kAz)  : fPointing->GetAz();
     433
     434    Double_t xsrcpos0   = fMatrix ? GetVal(kX)       : fSrcPos->GetX();
     435    Double_t ysrcpos0   = fMatrix ? GetVal(kY)       : fSrcPos->GetY();
     436    Double_t xmean0     = fMatrix ? GetVal(kMeanX)   : fHil->GetMeanX(); 
     437    Double_t ymean0     = fMatrix ? GetVal(kMeanY)   : fHil->GetMeanY();
     438    Double_t delta      = fMatrix ? GetVal(kDelta)   : fHil->GetDelta();
     439   
     440    Double_t size       = fMatrix ? GetVal(kSize)    : fHil->GetSize();
     441   
     442    Double_t m3long     = fMatrix ? GetVal(kM3Long)  : fHilExt->GetM3Long();
     443    Double_t asym       = fMatrix ? GetVal(kAsym)    : fHilExt->GetAsym();
    397444
    398445    //------------------------------------------
    399446    // convert to deg
    400     //    Double_t width    = width0  * fMm2Deg;
    401     //    Double_t length   = length0 * fMm2Deg;
    402447    Double_t xsrcpos  = xsrcpos0 * fMm2Deg;
    403448    Double_t ysrcpos  = ysrcpos0 * fMm2Deg;
    … …  
    453498    // Distance between estimated Disp and real source position
    454499    Double_t d2 = (xdisp-xsrcpos)*(xdisp-xsrcpos) +  (ydisp-ysrcpos)*(ydisp-ysrcpos);
    455     fHistChi2->Fill(d2);
     500    fHistMinPar->Fill(d2);
    456501   
    457502    // Longitudinal and transversal distances between Disp and real source positon
    … …  
    465510    fHistcosZA->Fill(cos(theta/kRad2Deg));
    466511   
    467     // to check the size dependence of the optimization
    468     fHistChi2Energy->Fill(log10(energy),d2);
    469     fHistChi2Size->Fill(log10(size),d2);
     512    // to check the size and energy dependence of the optimization
     513    fHistMinParEnergy->Fill(log10(energy),d2);
     514    fHistMinParSize->Fill(log10(size),d2);
    470515    fHistDuEnergy->Fill(log10(energy),Du);
    471516    fHistDuSize->Fill(log10(size),Du);
    … …  
    473518    fHistDvSize->Fill(log10(size),Dv);
    474519   
    475     // variables for Chi^2 calculation (= d^2)
     520    // variables for the Minimization parameter calculation (= d^2)
    476521    fNumEv += 1;
    477     fSumChi2 += d2; 
     522    fSumMinPar += d2; 
    478523   
    479524    return kTRUE;
    … …  
    483528// --------------------------------------------------------------------------
    484529//
    485 // Calculates the final Chi^2 of the Disp optimization
     530// Calculates the final Minimization parameter of the Disp optimization
    486531//
    487532Bool_t MHDisp::Finalize()
    488533{
    489     fChi2 = fSumChi2/fNumEv;
     534    fMinPar->SetVal(fSumMinPar/fNumEv);
    490535    *fLog << endl;
    491     *fLog << "MHDisp::Finalize: SumChi2, NumEv = " << fSumChi2 << ", " << fNumEv << endl;
    492     *fLog << "MHDisp::Finalize: Final Chi2 = " << fChi2 << endl;
    493    
    494     fChi2 = fHistChi2->GetMean();
    495     *fLog << "MHDisp::Finalize: Final Chi2 = " << fChi2 << endl;
     536    *fLog << "MHDisp::Finalize: SumMinPar, NumEv = " << fSumMinPar << ", " << fNumEv << endl;
     537    *fLog << "MHDisp::Finalize: Final MinPar = " << fMinPar->GetVal() << endl;
     538   
     539    fMinPar->SetVal(fHistMinPar->GetMean());
     540    *fLog << "MHDisp::Finalize: Final MinPar = " << fMinPar->GetVal() << endl;
    496541
    497542    SetReadyToSave();
    … …  
    542587    pad->cd(4);
    543588    gPad->SetBorderMode(0);
    544     fHistChi2->SetTitleOffset(1.2,"Y");
    545     fHistChi2->DrawClone(opt);
    546     fHistChi2->SetBit(kCanDelete);
    547     gPad->Modified();
    548 
    549     TProfile *profChi2Energy = fHistChi2Energy->ProfileX("Chi^2 vs. Energy",0,99999,"s");
    550     profChi2Energy->SetXTitle("log10(Energy) [GeV]");
    551     profChi2Energy->SetYTitle("Chi^2 = distance^2 Disp to real srcpos [deg^2]");
     589    fHistMinPar->SetTitleOffset(1.2,"Y");
     590    fHistMinPar->DrawClone(opt);
     591    fHistMinPar->SetBit(kCanDelete);
     592    gPad->Modified();
     593
     594    TProfile *profMinParEnergy = fHistMinParEnergy->ProfileX("Minimization parameter vs. Energy",0,99999,"s");
     595    profMinParEnergy->SetXTitle("log10(Energy) [GeV]");
     596    profMinParEnergy->SetYTitle("Minimization parameter = d^2 Disp to real srcpos [deg^2]");
    552597    pad->cd(5);
    553598    gPad->SetBorderMode(0);
    554     profChi2Energy->SetTitleOffset(1.2,"Y");
    555     profChi2Energy->DrawClone(opt);
    556     profChi2Energy->SetBit(kCanDelete);
    557     gPad->Modified();
    558 
    559     TProfile *profChi2Size = fHistChi2Size->ProfileX("Chi^2 vs. Size",0,99999,"s");
    560     profChi2Size->SetXTitle("log10(Size) [#phot]");
    561     profChi2Size->SetYTitle("Chi^2 = distance^2 Disp to real srcpos [deg^2]");
     599    profMinParEnergy->SetTitleOffset(1.2,"Y");
     600    profMinParEnergy->SetStats(0);
     601    profMinParEnergy->DrawClone(opt);
     602    profMinParEnergy->SetBit(kCanDelete);
     603    gPad->Modified();
     604
     605    TProfile *profMinParSize = fHistMinParSize->ProfileX("Minimization parameter vs. Size",0,99999,"s");
     606    profMinParSize->SetXTitle("log10(Size) [#phot]");
     607    profMinParSize->SetYTitle("Minimization parameter = d^2 Disp to real srcpos [deg^2]");
    562608    pad->cd(6);
    563609    gPad->SetBorderMode(0);
    564     profChi2Size->SetTitleOffset(1.2,"Y");
    565     profChi2Size->DrawClone(opt);
    566     profChi2Size->SetBit(kCanDelete);
     610    profMinParSize->SetTitleOffset(1.2,"Y");
     611    profMinParSize->SetStats(0);
     612    profMinParSize->DrawClone(opt);
     613    profMinParSize->SetBit(kCanDelete);
    567614    gPad->Modified();
    568615
    … …  
    604651    gPad->SetBorderMode(0);
    605652    profDuEnergy->SetTitleOffset(1.2,"Y");
     653    profDuEnergy->SetStats(0);
    606654    profDuEnergy->DrawClone(opt);
    607655    profDuEnergy->SetBit(kCanDelete);
    … …  
    614662    gPad->SetBorderMode(0);
    615663    profDuSize->SetTitleOffset(1.2,"Y");
     664    profDuSize->SetStats(0);
    616665    profDuSize->DrawClone(opt);
    617666    profDuSize->SetBit(kCanDelete);
    … …  
    634683    gPad->SetBorderMode(0);
    635684    profDvEnergy->SetTitleOffset(1.2,"Y");
     685    profDvEnergy->SetStats(0);
    636686    profDvEnergy->DrawClone(opt);
    637687    profDvEnergy->SetBit(kCanDelete);
    … …  
    644694    gPad->SetBorderMode(0);
    645695    profDvSize->SetTitleOffset(1.2,"Y");
     696    profDvSize->SetStats(0);
    646697    profDvSize->DrawClone(opt);
    647698    profDvSize->SetBit(kCanDelete);

  • trunk/MagicSoft/Mars/mtemp/mifae/library/MHDisp.h

    r5671 r5904  
    2121class MPointingPos;
    2222class MHMatrix;
     23class MParameterD;
    2324
    2425class MHDisp : public MH
    … …  
    3637    TString fImageParDispName;
    3738
    38     Int_t fSelectedPos;      // flag to select which Disp estimated position
    39                              // we want to fill the histograms (see Set function)
     39    Int_t fSelectedPos;      // flag to select which of the two Disp source position
     40                             // solutions we want to fill the histograms (see Set function)
    4041
    4142    TH1F *fHistEnergy;       // Energy distribution of events
    … …  
    4344    TH1F *fHistcosZA;        // cosinus Zenith angle distribution of events
    4445    TH2F *fSkymapXY;         // 2D histogram for Disp estimated source positions
    45     TH1F *fHistChi2;         // Histogram of the event Chi^2 (= d^2 = distance^2
    46                              // between Disp estimated and true source position)
     46    TH1F *fHistMinPar;       // Histogram of the event Minimization Parameter (= d^2 =
     47                             // = distance^2 between Disp estimated and true source position)
    4748    TH2F *fHistDuDv;         // Distribution of longitudinal (Du) and transversal
    4849                             // (Dv) distances between Disp and true source position   
    49     TH2F *fHistChi2Energy;   // Chi^2 (= d^2) vs. Energy
    50     TH2F *fHistChi2Size;     // Chi^2 (= d^2) vs. Size
     50    TH2F *fHistMinParEnergy; // Minimization Parameter (= d^2) vs. Energy
     51    TH2F *fHistMinParSize;   // Minimization Parameter (= d^2) vs. Size
    5152    TH2F *fHistDuEnergy;     // Du vs. Energy
    5253    TH2F *fHistDuSize;       // Du vs. Size
    … …  
    5859    Double_t fMm2Deg;        // conversion factor from mm to deg
    5960
    60     MHMatrix *fMatrix;       // matrix storing variables needed for Disp optimization
    61     TArrayI fMap;            // array storing the matrix mapping of columns defined
    62                              // in MDispCalc::InitMapping
     61    MHMatrix *fMatrix;       // matrix storing variables needed for the Disp optimization
     62    TArrayI fMap;            // array storing the matrix mapping column numbers corresponding
     63                             // to each variable added to fMatrix
    6364
    6465    Double_t GetVal(Int_t i) const;
    6566
    6667    Int_t fNumEv;            // total number of events
    67     Double_t fSumChi2;       // current sum of Chi^2      
    68     Double_t fChi2;          // final Chi^2 of the Disp optimization
     68    Double_t fSumMinPar;     // current sum of the minimization parameter   
     69    MParameterD *fMinPar;    // final minimization parameters of the Disp optimization
    6970
    7071public:
    … …  
    8182
    8283    void SetMatrixMap(MHMatrix *matrix, TArrayI &map);   
     84    void GetMatrixMap(MHMatrix* &matrix, TArrayI &map);  // get matrix and its mapping array
    8385
    84     Double_t GetChi2()     { return fChi2; }
     86    void InitMapping(MHMatrix *mat);       // define the matrix of variables
     87                                           // needed for the Disp optimization
    8588
    8689    TH1F *GetHistEnergy()       { return fHistEnergy; }
    … …  
    8891    TH1F *GetHistcosZA()        { return fHistcosZA; }
    8992    TH2F *GetSkymapXY()         { return fSkymapXY; }
    90     TH1F *GetHistChi2()         { return fHistChi2; }
     93    TH1F *GetHistMinPar()       { return fHistMinPar; }
    9194    TH2F *GetHistDuDv()         { return fHistDuDv; }
    92     TH2F *GetHistChi2Energy()   { return fHistChi2Energy; }
    93     TH2F *GetHistChi2Size()     { return fHistChi2Size; }
     95    TH2F *GetHistMinParEnergy() { return fHistMinParEnergy; }
     96    TH2F *GetHistMinParSize()   { return fHistMinParSize; }
    9497    TH2F *GetHistDuEnergy()     { return fHistDuEnergy; }
    9598    TH2F *GetHistDuSize()       { return fHistDuSize; }
    … …  
    100103    void Draw(Option_t *opt="");
    101104
    102     ClassDef(MHDisp, 1) // Container holding the Histograms for Disp
     105    ClassDef(MHDisp, 1) // Container holding the Histograms for Disp and
     106                        // the parameter to minimize in the Disp optimization
    103107};
    104108

  • trunk/MagicSoft/Mars/mtemp/mifae/library/Makefile

    r5672 r5904  
    7474        MTopologyCalc.cc \
    7575        MImageParDisp.cc \
    76         MDisp.cc \
     76        MDispParameters.cc \
    7777        MDispCalc.cc \
    7878        MHDisp.cc \
Note: See TracChangeset for help on using the changeset viewer.