Changeset 8656 for trunk/MagicSoft


Ignore:
Timestamp:
08/03/07 13:00:52 (17 years ago)
Author:
tbretz
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • trunk/MagicSoft/Mars/macros/optim/optimdisp.C

    r7409 r8656  
    88
    99    // -------------------- Setup ----------------------------
    10     opt.AddParameter("1-(MHillas.fWidth/MHillas.fLength)"); // M[0]
    11     opt.AddParameter("log10(MNewImagePar.fLeakage1+1)");    // M[1]
    1210
    13     opt.SetParameter(0,  1.30871);                          // Setup [0]
    14     opt.SetParameter(1,  5.81119);                          // Setup [1]
    15     opt.SetParameter(2,  0.763486);                         // Setup [2]
     11    opt.AddParameter("1-MHillas.fWidth/MHillas.fLength");   // M[0]
     12    opt.AddParameter("MHillasExt.fSlopeLong*sign(MHillasSrc.fCosDeltaAlpha)/MGeomCam.fConvMm2Deg"); // M[1]
     13    opt.AddParameter("MNewImagePar.fLeakage1");             // M[2]
     14    opt.AddParameter("log10(MHillas.fSize)");               // M[3]
    1615
    17     char *r = "([0]+([1]*pow(M[1], [2])))*M[0]";            // Rule to calc Disp
     16    // -------------- Parametrization --------------------
     17
     18    char *r = "M[0]*([0] + [1]*M[1] + [2]*M[2] + (M[3]>[3])*[4]*(M[3]-[3])^2)";
     19
     20    opt.FixParameter(0, 1.266195);
     21    opt.FixParameter(1, 0.100577);
     22    opt.FixParameter(2, 1.80309);
     23    opt.FixParameter(3, 2.87177);
     24    opt.FixParameter(4, 0.616823);
     25
     26    opt.AddPreCut("MNewImagePar.fLeakage1>0");
     27    opt.AddPreCut("log10(MHillas.fSize)<3.2");
     28
     29    //opt.AddPreCut("abs(MHillasSrc.fDCA*MGeomCam.fConvMm2Deg)<0.2");
     30    //opt.AddPreCut("(MHillasSrc.fDist*MGeomCam.fConvMm2Deg-0.5)*7.2>MHillasExt.fSlopeLong*sign(MHillasSrc.fCosDeltaAlpha)/MGeomCam.fConvMm2Deg");
     31    opt.AddPreCut("MHillasExt.fM3Long*sign(MHillasSrc.fCosDeltaAlpha)*MGeomCam.fConvMm2Deg>-0.07");
     32    opt.AddPreCut("DataType.fVal>0.5");
    1833
    1934    // -------------------- Run ----------------------------
    … …  
    2641     MFMagicCuts cuts;
    2742     cuts.SetHadronnessCut(MFMagicCuts::kArea);
    28      cuts.SetThetaCut(MFMagicCuts::kOn);
     43     cuts.SetThetaCut(MFMagicCuts::kNone);
    2944
    3045     TArrayD arr(10);
    3146     arr[0]=  1.3245;
    3247     arr[1]=  0.208700;
    33      arr[2]=  0.229200;
    34      arr[3]=  5.305200;
    35      arr[4]=  0.098930;
    36      arr[5]= -0.082950;
    37      arr[6]=  8.2957;
    38      arr[7]=  0.8677;
     48     arr[2]=  0.0836169;
     49     arr[3]=  5.63973;
     50     arr[4]=  0.0844195;
     51     arr[5]= -0.07;
     52     arr[6]=  13.2;
     53     arr[7]=  1.0;
     54     arr[8]=  7.2;
     55     arr[9]=  0.5;
    3956
    4057     cuts.SetVariables(arr);
    … …  
    4764
    4865     -------------------- Other cuts ----------------------
    49      opt.AddPreCut("MNewImagePar.fLeakage1<0.0001");
    50     */
     66     */
    5167
    52     opt.RunDisp("ganymedmcpart.root", r);
     68    opt.RunDisp("ganymed00000001-summary.root", r);
    5369}
    5470
    5571/* ------------------ Good strategy -------------------
    56  1) first fix parameters:
    57       opt.FixParameter(1,  0);
    58       opt.FixParameter(2,  1);
    59     and process only showers without leakage
    60       opt.AddPreCut("MNewImagePar.fLeakage1<0.0001");
    61  2) release parameters 1 and 2 and fix 0 to the result of 1)
    62       opt.FixParameter(0,  0.8362);
    63       opt.SetParameter(1,  2.0);
    64       opt.SetParameter(2,  0.8);
    65     and process only showers with leakage
    66       opt.AddPreCut("MNewImagePar.fLeakage1>0.0001");
    67  3) release all parameters and start with the result of 1) and 2)
    68       opt.SetParameter(0,  0.8362);
    69       opt.SetParameter(1,  5.84);
    70       opt.SetParameter(2,  0.76);
    71     and process all showers. (Comment out opt.PreCuts())
    72  */
     72
     73   Par  |   0    1     2      3     4     |  Cut
     74 -------+---------------------------------+-----------------------
     75  Fit 1 |  1.3  0.8  fix=0  fix=0  fix=0  | Leak1==0 lgSize<2.5
     76  Fit 2 |  fix  fix  fix=0   2.4    0.3   | Leak1==0
     77  Fit 2 |  fix  fix   1.8    fix    fix   | Leak1>0
     78  Fit 3 |  free free  fix    fix    fix   | -/-
     79 -------+---------------------------------+-----------------------
     80*/

  • trunk/MagicSoft/Mars/mjtrain/MJTrainDisp.cc

    r8644 r8656  
    189189// Run Disp optimization
    190190//
    191 Bool_t MJTrainDisp::Train(const char *out, const MDataSet &set, Int_t num)
     191Bool_t MJTrainDisp::TrainDisp(const char *out, const MDataSet &set, Int_t num)
    192192{
    193193    SetTitle(Form("TrainDisp: %s", out));
    … …  
    224224    if (fEnableWeights)
    225225        train.AddColumn("MWeight.fVal");
    226     train.AddColumn(fgTrainParameter);
     226    train.AddColumn(fTrainParameter);
    227227
    228228    // ----------------------- Fill Matrix RF ----------------------
    … …  
    240240
    241241    // ------------------------ Train RF --------------------------
    242     MRanForestCalc rf(fTitle);
     242    MRanForestCalc rf("TrainDisp", fTitle);
    243243    rf.SetNumTrees(fNumTrees);
    244244    rf.SetNdSize(fNdSize);
    … …  
    268268    gLog.Separator("Test");
    269269
     270    MH::SetPalette("pretty");
     271
    270272    MParList  plist;
    271273    MTaskList tlist;
    … …  
    329331    MFillH fillh2a(&hdisp1, "", "FillSize");
    330332    MFillH fillh2b(&hdisp2, "", "FillEnergy");
    331     fillh2a.SetDrawOption("blue profx");
    332     fillh2b.SetDrawOption("blue profx");
     333    fillh2a.SetDrawOption("colz profx");
     334    fillh2b.SetDrawOption("colz profx");
    333335    fillh2a.SetNameTab("Size");
    334336    fillh2b.SetNameTab("Energy");
    … …  
    354356
    355357    // Print the result
    356     *fLog << inf << "Rule: " << rule << endl;
     358    *fLog << inf;
     359    *fLog << "Rule: " << rule << endl;
     360    *fLog << "Disp: " << fTrainParameter << endl;
    357361    hist.GetAlphaFitter().Print("result");
    358362
    … …  
    366370}
    367371
     372/*
     373#include "MParameterCalc.h"
     374#include "MHillasCalc.h"
     375#include "../mpointing/MSrcPosRndm.h"
     376
     377Bool_t MJTrainDisp::TrainGhostbuster(const char *out, const MDataSet &set, Int_t num)
     378{
     379    SetTitle(Form("TrainGhostbuster: %s", out));
     380
     381    if (fDisplay)
     382        fDisplay->SetTitle(fTitle);
     383
     384    if (!set.IsValid())
     385    {
     386        *fLog << err << "ERROR - DataSet invalid!" << endl;
     387        return kFALSE;
     388    }
     389
     390    if (!HasWritePermission(out))
     391        return kFALSE;
     392
     393    *fLog << inf;
     394    fLog->Separator(GetDescriptor());
     395
     396    // --------------------- Setup files --------------------
     397    MReadMarsFile readtrn("Events");
     398    MReadMarsFile readtst("Events");
     399    readtrn.DisableAutoScheme();
     400    readtst.DisableAutoScheme();
     401
     402    if (!set.AddFilesOn(readtrn))
     403        return kFALSE;
     404    if (!set.AddFilesOff(readtst))
     405        return kFALSE;
     406
     407    // ----------------------- Setup Matrix ------------------
     408    MHMatrix train("Train");
     409    train.AddColumns(fRules);
     410    if (fEnableWeights)
     411        train.AddColumn("MWeight.fVal");
     412    train.AddColumn("sign(MHillasSrc.fCosDeltaAlpha)==sign(SignStore.fVal)");
     413
     414    MParameterCalc calc("MHillasSrc.fCosDeltaAlpha", "SignStore");
     415    calc.SetNameParameter("SignStore");
     416
     417    MSrcPosRndm rndm;
     418    rndm.SetRule(fTrainParameter);
     419    //rndm.SetDistOfSource(120*3.37e-3);
     420
     421    MHillasCalc hcalc;
     422    hcalc.SetFlags(MHillasCalc::kCalcHillasSrc);
     423
     424    // ----------------------- Fill Matrix RF ----------------------
     425    MTFillMatrix fill(fTitle);
     426    fill.SetDisplay(fDisplay);
     427    fill.SetLogStream(fLog);
     428    fill.SetDestMatrix1(&train, num);
     429    fill.SetReader(&readtrn);
     430    fill.AddPreCuts(fPreCuts);
     431    fill.AddPreCuts(fTrainCuts);
     432    fill.AddPreTasks(fPreTasks);
     433    fill.AddPostTasks(fPostTasks);
     434    fill.AddPostTask(&calc);
     435    fill.AddPostTask(&rndm);
     436    fill.AddPostTask(&hcalc);
     437    if (!fill.Process())
     438        return kFALSE;
     439
     440    // ------------------------ Train RF --------------------------
     441    MRanForestCalc rf("TrainGhostbuster", fTitle);
     442    rf.SetNumTrees(fNumTrees);
     443    rf.SetNdSize(fNdSize);
     444    rf.SetNumTry(fNumTry);
     445    rf.SetNumObsoleteVariables(1);
     446    rf.SetLastDataColumnHasWeights(fEnableWeights);
     447    rf.SetDisplay(fDisplay);
     448    rf.SetLogStream(fLog);
     449    rf.SetFileName(out);
     450    rf.SetDebug(fDebug>1);
     451    rf.SetNameOutput("Sign");
     452
     453    if (!rf.TrainRegression(train))                  // regression (best choice)
     454        return kFALSE;
     455
     456    // --------------------- Display result ----------------------
     457
     458    gLog.Separator("Test");
     459
     460    MH::SetPalette("pretty");
     461
     462    MParList  plist;
     463    MTaskList tlist;
     464    plist.AddToList(this);
     465    plist.AddToList(&tlist);
     466    //plist.AddToList(&b);
     467
     468    MFilterList list;
     469    if (!list.AddToList(fPreCuts))
     470        *fLog << err << "ERROR - Calling MFilterList::AddToList for fPreCuts failed!" << endl;
     471    if (!list.AddToList(fTestCuts))
     472        *fLog << err << "ERROR - Calling MFilterList::AddToList for fTestCuts failed!" << endl;
     473
     474    MContinue cont(&list);
     475    cont.SetInverted();
     476
     477    const char *rule = "abs(Sign.fVal-(sign(MHillasSrc.fCosDeltaAlpha)==sign(SignStore.fVal)))";
     478
     479    //MChisqEval eval;
     480    //eval.SetY1("sqrt(ThetaSquared.fVal)");
     481
     482    MH3 hsign1("MHillas.fSize",  rule);
     483    MH3 hsign2("MMcEvt.fEnergy", rule);
     484    hsign1.SetTitle("Was ist das? vs. Size:Size [phe]:XXX");
     485    hsign2.SetTitle("Was ist das? vs. Energy:Enerhy [GeV]:XXXX");
     486
     487    MBinning binsx( 70, 10, 31623, "BinningMH3X", "log");
     488    MBinning binsy( 51,  0,     1, "BinningMH3Y", "lin");
     489
     490    plist.AddToList(&binsx);
     491    plist.AddToList(&binsy);
     492
     493    MFillH fillh2a(&hsign1, "", "FillSize");
     494    MFillH fillh2b(&hsign2, "", "FillEnergy");
     495    fillh2a.SetDrawOption("profx colz");
     496    fillh2b.SetDrawOption("profx colz");
     497    fillh2a.SetNameTab("Size");
     498    fillh2b.SetNameTab("Energy");
     499
     500    tlist.AddToList(&readtst);
     501    tlist.AddToList(&cont);
     502    tlist.AddToList(&calc);
     503    tlist.AddToList(&rndm);
     504    tlist.AddToList(&hcalc);
     505    tlist.AddToList(&rf);
     506    tlist.AddToList(&fillh2a);
     507    tlist.AddToList(&fillh2b);
     508    //tlist.AddToList(&eval);
     509
     510    MEvtLoop loop(fTitle);
     511    loop.SetLogStream(fLog);
     512    loop.SetDisplay(fDisplay);
     513    loop.SetParList(&plist);
     514    //if (!SetupEnv(loop))
     515    //    return kFALSE;
     516
     517    if (!loop.Eventloop())
     518        return kFALSE;
     519
     520    // Print the result
     521    *fLog << inf << "Rule: " << rule << endl;
     522
     523    //DisplayResult(hdisp1, hdisp2);
     524
     525    SetPathOut(out);
     526    if (!WriteDisplay(0, "UPDATE"))
     527        return kFALSE;
     528
     529    return kTRUE;
     530}
     531*/

  • trunk/MagicSoft/Mars/mjtrain/MJTrainDisp.h

    r8644 r8656  
    2727
    2828    Bool_t Train(const char *out, const MDataSet &set, Int_t num);
     29    //Bool_t TrainGhostbuster(const char *out, const MDataSet &set, Int_t num);
    2930
    3031    ClassDef(MJTrainDisp, 0)//Class to train Random Forest disp estimator
Note: See TracChangeset for help on using the changeset viewer.