/* More explanations at the end of the file */ void optimdisp() { MJOptimizeDisp opt; // Initialize optimization class opt.SetDebug(2); // Choose the level of output opt.SetOptimizer(MJOptimize::kMigrad); // Choose the fit algorithm (see MJOptimize) opt.EnableTestTrain(); // Split sample into test and train // -------------------- Setup ---------------------------- opt.AddParameter("1-(MHillas.fWidth/MHillas.fLength)"); // M[0] opt.AddParameter("log10(MNewImagePar.fLeakage1+1)"); // M[1] opt.SetParameter(0, 1.30871); // Setup [0] opt.SetParameter(1, 5.81119); // Setup [1] opt.SetParameter(2, 0.763486); // Setup [2] char *r = "([0]+([1]*pow(M[1], [2])))*M[0]"; // Rule to calc Disp // -------------------- Run ---------------------------- MStatusDisplay *d = new MStatusDisplay; opt.SetDisplay(d); /* -------------------- Magic-Cuts ---------------------- MFMagicCuts cuts; cuts.SetHadronnessCut(MFMagicCuts::kArea); cuts.SetThetaCut(MFMagicCuts::kOn); TArrayD arr(10); arr[0]= 1.3245; arr[1]= 0.208700; arr[2]= 0.229200; arr[3]= 5.305200; arr[4]= 0.098930; arr[5]= -0.082950; arr[6]= 8.2957; arr[7]= 0.8677; cuts.SetVariables(arr); opt.AddPreCut(&cuts); -------------------- Energy Slope -------------------- MFEnergySlope slope(-2.8); opt.AddPreCut(&slope); -------------------- Other cuts ---------------------- opt.AddPreCut("MNewImagePar.fLeakage1<0.0001"); */ opt.RunDisp("ganymedmcpart.root", r); } /* ------------------ Good strategy ------------------- 1) first fix parameters: opt.FixParameter(1, 0); opt.FixParameter(2, 1); and process only showers without leakage opt.AddPreCut("MNewImagePar.fLeakage1<0.0001"); 2) release parameters 1 and 2 and fix 0 to the result of 1) opt.FixParameter(0, 0.8362); opt.SetParameter(1, 2.0); opt.SetParameter(2, 0.8); and process only showers with leakage opt.AddPreCut("MNewImagePar.fLeakage1>0.0001"); 3) release all parameters and start with the result of 1) and 2) opt.SetParameter(0, 0.8362); opt.SetParameter(1, 5.84); opt.SetParameter(2, 0.76); and process all showers. (Comment out opt.PreCuts()) */