Index: trunk/FACT++/spectrum/display.C =================================================================== --- trunk/FACT++/spectrum/display.C (revision 19970) +++ trunk/FACT++/spectrum/display.C (revision 19971) @@ -91,4 +91,21 @@ h1->SetMarkerStyle(23); h1->DrawCopy("P same"); + + // -------------------------------- + + c = new TCanvas("Impact", "Impact"); + c->Divide(2,2); + c->cd(1); + file.GetObject("MC/dense/TrueEnergy/Impact", h2); + h2->DrawCopy("colz"); + h2->DrawCopy("same"); + c->cd(2); + file.GetObject("MC/sparse/TrueEnergy/Impact", h2); + h2->DrawCopy("colz"); + h2->DrawCopy("same"); + c->cd(3); + file.GetObject("MC/theta/Impact", h2); + h2->DrawCopy("colz"); + h2->DrawCopy("same"); // -------------------------------- Index: trunk/FACT++/spectrum/simulation.sql =================================================================== --- trunk/FACT++/spectrum/simulation.sql (revision 19970) +++ trunk/FACT++/spectrum/simulation.sql (revision 19971) @@ -10,6 +10,7 @@ INTERVAL(LOG10(Energy), %108:sparse) AS `.sparse_sim`, INTERVAL(LOG10(Energy), %109:dense) AS `.dense_sim`, + INTERVAL(Impact/100, %110:impact) AS `.impact`, - (%110:spectrum)/POW(Energy, SpectralIndex) AS SpectralWeight, -- FIXME: Is this correct for files with different Slopes? + (%111:spectrum)/POW(Energy, SpectralIndex) AS SpectralWeight, -- FIXME: Is this correct for files with different Slopes? LogEnergyEst - log10(Energy) AS Residual FROM @@ -27,4 +28,5 @@ `.dense_est`, `.dense_sim`, + `.impact`, -- Without any weight applied @@ -59,3 +61,3 @@ ThetaDist USING(`.theta`) GROUP BY - `.theta`, `.sparse_est`, `.sparse_sim`, `.dense_est`, `.dense_sim` + `.theta`, `.sparse_est`, `.sparse_sim`, `.dense_est`, `.dense_sim`, `.impact` Index: trunk/FACT++/src/spectrum.cc =================================================================== --- trunk/FACT++/src/spectrum.cc (revision 19970) +++ trunk/FACT++/src/spectrum.cc (revision 19971) @@ -139,10 +139,12 @@ po::options_description binnings("Binnings"); binnings.add_options() - ("theta", var(Binning(90, 0, 90)), "Add equidistant bins in theta (degrees). Syntax: N,lo,hi (Number N of equidistant bins between lo and hi)") - ("theta-bin", vars(), "Add a bin-edge to the theta binning (degree)") - ("energy-dense", var(Binning(30, 2, 5)), "Add equidistant bins in log10 simulated energy. Syntax: N,lo,hi (Number N of equidistant bins between lo and hi)") - ("energy-dense-bin", vars(), "Add a bin-edge to the binnnig in log10 simulated enegry") - ("energy-sparse", var(Binning(15, 2, 5)), "Add equidistant bins in log10 estimated energy. Syntax: N,lo,hi (Number N of equidistant bins between lo and hi)") - ("energy-sparse-bin", vars(), "Add a bin-edge to the binning in log10 estimated enegry") + ("theta", var(Binning(90, 0, 90)), "Add equidistant bins in theta (degrees). Syntax: N,lo,hi (Number N of equidistant bins between lo and hi)") + ("theta-bin", vars(), "Add a bin-edge to the theta binning (degree)") + ("energy-dense", var(Binning(30, 2, 5)), "Add equidistant bins in log10 simulated energy. Syntax: N,lo,hi (Number N of equidistant bins between lo and hi)") + ("energy-dense-bin", vars(), "Add a bin-edge to the binnnig in log10 simulated enegry") + ("energy-sparse", var(Binning(15, 2, 5)), "Add equidistant bins in log10 estimated energy. Syntax: N,lo,hi (Number N of equidistant bins between lo and hi)") + ("energy-sparse-bin", vars(), "Add a bin-edge to the binning in log10 estimated enegry") + ("impact", var(Binning(28, 0, 280)), "Add equidistant bins in impact in meter. Syntax: N,lo,hi (Number N of equidistant bins between lo and hi)") + ("impact-bin", vars(), "Add a bin-edge to the binning in impact in meter") ; @@ -632,9 +634,11 @@ Binning binning_dense = conf.Get("energy-dense"); Binning binning_sparse = conf.Get("energy-sparse"); + Binning binning_impact = conf.Get("impact"); cout << '\n'; - cout << "Binning 'theta': " << binning_theta.str() << endl; + cout << "Binning 'theta': " << binning_theta.str() << endl; cout << "Binning 'dense': " << binning_dense.str() << endl; - cout << "Binning 'sparse': " << binning_sparse.str() << endl; + cout << "Binning 'sparse': " << binning_sparse.str() << endl; + cout << "Binning 'impact': " << binning_impact.str() << endl; const uint16_t source_key = conf.Get("source-key"); @@ -656,4 +660,5 @@ const string str_dense = binning_dense.list(); const string str_sparse = binning_sparse.list(); + const string str_impact = binning_impact.list(); // ------------------------------------------------------------------------- @@ -745,8 +750,10 @@ CreateBinning(connection, qlog, binning_dense, "Energy_dense", "Dense binning in log10 Energy"); CreateBinning(connection, qlog, binning_sparse, "Energy_sparse", "Sparse binning in log10 Energy"); + CreateBinning(connection, qlog, binning_impact, "Impact", "Binning in impact distance"); Dump(flog, connection, "BinningTheta"); Dump(flog, connection, "BinningEnergy_dense"); Dump(flog, connection, "BinningEnergy_sparse"); + Dump(flog, connection, "BinningImpact"); // ------------------------------------------------------------------- @@ -928,6 +935,6 @@ " LEFT JOIN\n" " factmc.CorsikaSetup ON FileId=RUNNR\n" - " GROUP BY\n" - " `CSCAT[1]`, `CSCAT[2]`\n" + // " GROUP BY\n" + // " `CSCAT[1]`, `CSCAT[2]`\n" " ORDER BY\n" " MaxImpactHi\n" @@ -955,5 +962,4 @@ cout << "Execution time: " << sec4 << "s\n\n"; } - // ------------------------------------------------------------------- @@ -1297,5 +1303,6 @@ " INDEX (`.dense_sim`) USING HASH,\n" " INDEX (`.sparse_est`, `.sparse_sim`) USING HASH,\n" - " INDEX (`.dense_est`, `.dense_sim`) USING HASH\n" + " INDEX (`.dense_est`, `.dense_sim`) USING HASH,\n" + " INDEX (`.impact`) USING HASH\n" ") ENGINE=MEMORY COMMENT='Sum of counts and (squared) weightes of Monte Carlo Data after analysis'\n" "AS\n" @@ -1317,5 +1324,5 @@ //query6.template_defaults["columns"] = "FileId, EvtNumber, CorsikaNumReuse,"; - query8.template_defaults["columns"] = "Energy, SpectralIndex,"; + query8.template_defaults["columns"] = "Energy, SpectralIndex, Impact,"; query8.template_defaults["zenith"] = "DEGREES(Theta)"; query8.template_defaults["files"] = "MonteCarloFiles"; @@ -1327,4 +1334,5 @@ query8.template_defaults["dense"] = str_dense.c_str(); query8.template_defaults["theta"] = str_theta.c_str(); + query8.template_defaults["impact"] = str_impact.c_str(); query8.template_defaults["spectrum"] = spectrum.c_str(); query8.template_defaults["estimator"] = estimator.c_str(); @@ -1619,6 +1627,72 @@ WriteHistogram(connection, hist_sim); + // ------------------------------------------------------------------- + // ------------------------ ImpactDistribution ----------------------- + // ------------------------------------------------------------------- + + cout << separator("ImpactDistribution_"+*ib) << '\n'; + + Time start13; + + mysqlpp::Query query13(&connection); + query13 << + "CREATE TEMPORARY TABLE ImpactDistribution_%100:binning\n" + "(\n" + " SignalN INT UNSIGNED NOT NULL\n" + ") ENGINE=MEMORY COMMENT='Impact Distribution of unweighted signal events after cuts'\n" + "AS\n" + "(\n" + " SELECT\n" + " %101:bin,\n" + " `.impact`,\n" + " SUM(SignalN) AS SignalN\n" + " FROM\n" + " AnalysisMC\n" + " GROUP BY\n" + " %101:bin, `.impact`\n" + " ORDER BY\n" + " %101:bin, `.impact`\n" + ")"; + + query13.parse(); + + query13.template_defaults["binning"] = ib->c_str(); + query13.template_defaults["bin"] = *ib=="theta" ? "`.theta`" : ("`."+*ib+"_sim`").c_str(); + + if (print_queries) + PrintQuery(query13.str()); + + qlog << query13 << ";\n" << endl; + if (connection.connected()) + { + cout << query13.execute().info() << endl; + ShowWarnings(connection); + Dump(flog, connection, "ImpactDistribution_"+*ib); + + const auto sec13 = Time().UnixTime()-start13.UnixTime(); + cout << "Execution time: " << sec13 << "s\n"; + } + + hist_sim.name = "Impact"; + hist_sim.title = "Distribution of Impact Distance"; + hist_sim.y = ".impact"; + hist_sim.v = "SignalN"; + hist_sim.err = ""; + hist_sim.table = "ImpactDistribution_"+*ib; + hist_sim.binningy = binning_impact; + hist_sim.axisy = "Impact Distance [m]"; + hist_sim.axisz = "Counts"; + + WriteHistogram(connection, hist_sim); + + // =================================================================== + // =================================================================== + if (*ib=="theta") continue; + + // =================================================================== + // =================================================================== + // ------------------------------------------------------------------- @@ -1796,5 +1870,5 @@ // ------------------------------------------------------------------- - // ------------------------- SimulatedSpectrum ----------------------- + // -------------------------- MigrationMatrix ------------------------ // -------------------------------------------------------------------