Changeset 19913 for trunk/FACT++


Ignore:
Timestamp:
12/15/19 23:30:18 (5 years ago)
Author:
tbretz
Message:
Write spectrum.C again.
File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/FACT++/src/spectrum.cc

    r19911 r19913  
    147147        ;
    148148
    149     po::options_description queries("Analysis Query");
    150     queries.add_options()
     149    po::options_description analysis("Analysis Setup");
     150    analysis.add_options()
    151151        ("analysis",    var<string>("analysis.sql"),   "File with the analysis query. A default file is created automatically in the <prefix> directory it does not exist.")
    152         //("data",        var<string>("data.sql"),       "File with the query which creates the data summary. A default file is created automatically in the <prefix> directory it does not exist.")
    153         //("simulation",  var<string>("simulation.sql"), "File with the query which creates the simulation summary. A default file is created automatically in the <prefix> directory it does not exist.")
    154         ;
    155 
    156     po::options_description preparation("Preparation");
    157     preparation.add_options()
    158         ("source-key", var<uint16_t>(5),          "Source key to be used in file selection.")
    159         ("selector",   vars<string>(),            "WHERE clause to be used in file selection.")
     152        ("source-key", var<uint16_t>(5),          "Source key to be used in data file selection.")
     153        ("selector",   vars<string>(),            "WHERE clause to be used in data file selection.")
    160154        ("estimator",  var<string>()->required(), "Energy estimator to be used.")
    161155        ("spectrum",   var<string>()->required(), "Spectral shape for re-weighting of simulated 'Energy'")
    … …  
    181175    conf.AddOptions(control);
    182176    conf.AddOptions(binnings);
    183     conf.AddOptions(queries);
    184     conf.AddOptions(preparation);
     177    conf.AddOptions(analysis);
    185178    conf.AddOptions(debug);
    186179    //conf.SetArgumentPositions(p);
    … …  
    369362
    370363    const auto vec = bins.vec();
    371     for (size_t i=1; i<vec.size()-2; i++)
     364    for (size_t i=1; i<vec.size()-1; i++)
    372365        query1 << "  ( " << i << ", " << vec[i-1] << ", " << vec[i] << " ),\n";
    373     query1 << "  ( " << vec.size()-2 << ", " << vec[vec.size()-2] << ", " << vec[vec.size()-1] << " )\n";
     366    query1 << "  ( " << vec.size()-1 << ", " << vec[vec.size()-2] << ", " << vec[vec.size()-1] << " )\n";
    374367
    375368    qlog << query1 << ";\n" << endl;
    … …  
    787780
    788781    query1.parse();
    789     for (auto it=env.cbegin(); it!=env.cend(); it++)
    790         query1.template_defaults[it->first.c_str()] = it->second.c_str();
     782    //for (auto it=env.cbegin(); it!=env.cend(); it++)
     783    //    query1.template_defaults[it->first.c_str()] = it->second.c_str();
    791784
    792785    query1.template_defaults["source"] = to_string(source_key).c_str();
    … …  
    842835
    843836    query2.parse();
    844     for (auto it=env.cbegin(); it!=env.cend(); it++)
    845         query2.template_defaults[it->first.c_str()] = it->second.c_str();
     837    //for (auto it=env.cbegin(); it!=env.cend(); it++)
     838    //    query2.template_defaults[it->first.c_str()] = it->second.c_str();
    846839
    847840    query2.template_defaults["bins"] = str_theta.c_str();
    … …  
    896889
    897890    query3.parse();
    898     for (auto it=env.cbegin(); it!=env.cend(); it++)
    899         query3.template_defaults[it->first.c_str()] = it->second.c_str();
     891    //for (auto it=env.cbegin(); it!=env.cend(); it++)
     892    //    query3.template_defaults[it->first.c_str()] = it->second.c_str();
    900893
    901894    if (print_queries)
    … …  
    942935
    943936    query4.parse();
    944     for (auto it=env.cbegin(); it!=env.cend(); it++)
    945         query4.template_defaults[it->first.c_str()] = it->second.c_str();
     937    //for (auto it=env.cbegin(); it!=env.cend(); it++)
     938    //    query4.template_defaults[it->first.c_str()] = it->second.c_str();
    946939
    947940    if (print_queries)
    … …  
    10171010
    10181011    query5.parse();
    1019     for (auto it=env.cbegin(); it!=env.cend(); it++)
    1020         query5.template_defaults[it->first.c_str()] = it->second.c_str();
     1012    //for (auto it=env.cbegin(); it!=env.cend(); it++)
     1013    //    query5.template_defaults[it->first.c_str()] = it->second.c_str();
    10211014
    10221015    query5.template_defaults["table"]    = "OriginalMC";
    … …  
    11111104
    11121105    query6.parse();
    1113     for (auto it=env.cbegin(); it!=env.cend(); it++)
    1114         query6.template_defaults[it->first.c_str()] = it->second.c_str();
     1106    //for (auto it=env.cbegin(); it!=env.cend(); it++)
     1107    //    query6.template_defaults[it->first.c_str()] = it->second.c_str();
    11151108
    11161109    if (print_queries)
    … …  
    12121205
    12131206    query7.parse();
    1214     for (auto it=env.cbegin(); it!=env.cend(); it++)
    1215         query7.template_defaults[it->first.c_str()] = it->second.c_str();
     1207    //for (auto it=env.cbegin(); it!=env.cend(); it++)
     1208    //    query7.template_defaults[it->first.c_str()] = it->second.c_str();
    12161209
    12171210    query7.template_defaults["table"] = "OriginalMC";
    … …  
    14821475
    14831476        query9.parse();
    1484         for (auto it=env.cbegin(); it!=env.cend(); it++)
    1485             query9.template_defaults[it->first.c_str()] = it->second.c_str();
     1477        //for (auto it=env.cbegin(); it!=env.cend(); it++)
     1478        //    query9.template_defaults[it->first.c_str()] = it->second.c_str();
    14861479
    14871480        query9.template_defaults["table"]    = table.c_str();
    … …  
    17131706
    17141707        query10.parse();
    1715         for (auto it=env.cbegin(); it!=env.cend(); it++)
    1716             query10.template_defaults[it->first.c_str()] = it->second.c_str();
     1708        //for (auto it=env.cbegin(); it!=env.cend(); it++)
     1709        //    query10.template_defaults[it->first.c_str()] = it->second.c_str();
    17171710
    17181711        query10.template_defaults["binning"] = ib->c_str();
    … …  
    18311824
    18321825        query11.parse();
    1833         for (auto it=env.cbegin(); it!=env.cend(); it++)
    1834             query11.template_defaults[it->first.c_str()] = it->second.c_str();
     1826        //for (auto it=env.cbegin(); it!=env.cend(); it++)
     1827        //    query11.template_defaults[it->first.c_str()] = it->second.c_str();
    18351828
    18361829        query11.template_defaults["binning"] = ib->c_str();
    … …  
    19431936    // --------------------------- Spectrum ------------------------------
    19441937    // -------------------------------------------------------------------
     1938
     1939    sindent mindent(mlog);
     1940    mlog << "void spectrum()\n";
     1941    mlog << "{\n" << indent(4);
     1942    mlog <<
     1943        "TGraphErrors *g = new TGraphErrors;\n"
     1944        "g->SetMarkerStyle(kFullDotMedium);\n\n"
     1945        "TGraph *ul = new TGraph;\n"
     1946        "ul->SetMarkerStyle(23);\n\n";
    19451947
    19461948    const vector<string> targets = { "Theta", "Energy" };
    … …  
    20012003
    20022004        query13.parse();
    2003         for (auto it=env.cbegin(); it!=env.cend(); it++)
    2004             query13.template_defaults[it->first.c_str()] = it->second.c_str();
     2005        //for (auto it=env.cbegin(); it!=env.cend(); it++)
     2006        //    query13.template_defaults[it->first.c_str()] = it->second.c_str();
    20052007
    20062008        query13.template_defaults["table"]     = table.c_str();
    … …  
    20742076                const size_t bin = row[*ib=="Theta" ? ".theta" : ".sparse_est"];
    20752077
     2078                const double flux   = row["Flux"];
     2079                const double error  = row["ErrFlux"];
     2080                const double center = row["center"];
     2081                const double sigma  = row["SigmaFlux"];
     2082
     2083                if (*ib=="Energy" && flux>0)
     2084                {
     2085                    mlog << "g->SetPoint(g->GetN(), pow(10, " << center << "), " << flux << ");\n";
     2086                    mlog << "g->SetPointError(g->GetN()-1, 0, " << error << ");\n";
     2087                }
     2088
    20762089#ifdef HAVE_ROOT
    20772090                const double dat_sig  = row["Signal"];
    … …  
    20982111                rolke.SetPoissonBkgKnownEff(dat_isig, dat_ibg*5, 5, ieff);
    20992112                rolke_int[bin] = rolke.GetUpperLimit()/areatime;
     2113
     2114                if (*ib=="Energy" && (sigma<1 || dat_sig<10 || dat_bg<2))
     2115                    mlog << "ul->SetPoint(ul->GetN(), pow(10, " << center << "), " << double(rolke_ul[bin]) << ");\n";
    21002116#endif
     2117
    21012118                if (verbose>0)
    21022119                {
    2103                     cout << setw(5)  << row["center"] << ":";
     2120                    cout << setw(5) << center << ":";
    21042121                    cout << " " << setw(10) << row["Excess"];
    21052122                    cout << " " << setw(10) << row["Significance"];
    2106                     cout << " " << setw(10) << row["Flux"];
    2107                     cout << " " << setw(10) << row["ErrFlux"];
     2123                    cout << " " << setw(10) << flux;
     2124                    cout << " " << setw(10) << error;
    21082125                    cout << endl;
    21092126                }
    … …  
    22322249    }
    22332250
     2251    mlog << "\n"
     2252        //"g.DrawClone(\"AP\");\n"
     2253        //"ul.DrawClone(\"P\");\n\n"
     2254        "TMultiGraph mg;\n"
     2255        "mg.SetTitle(\"Differential Energy Spectrum;E [GeV];dN/dE [cm^{-2} s^{-1} TeV^{-1}]\");\n"
     2256        "mg.Add(g,  \"P\");\n"
     2257        "mg.Add(ul, \"P\");\n"
     2258        "mg.DrawClone(\"A\");\n\n"
     2259        "gPad->SetLogx();\n"
     2260        "gPad->SetLogy();\n";
     2261    mlog << indent(0) << "}" << endl;
     2262
    22342263    // -------------------------------------------------------------------
    22352264    // ----------------------------- Summary -----------------------------
Note: See TracChangeset for help on using the changeset viewer.