Changeset 19971


Ignore:
Timestamp:
Jul 20, 2020, 4:20:59 PM (4 weeks ago)
Author:
tbretz
Message:
Added a binning in Impact Parameter for crosschecks
Location:
trunk/FACT++
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • trunk/FACT++/spectrum/display.C

    r19909 r19971  
    9191    h1->SetMarkerStyle(23);
    9292    h1->DrawCopy("P same");
     93
     94    // --------------------------------
     95
     96    c = new TCanvas("Impact", "Impact");
     97    c->Divide(2,2);
     98    c->cd(1);
     99    file.GetObject("MC/dense/TrueEnergy/Impact", h2);
     100    h2->DrawCopy("colz");
     101    h2->DrawCopy("same");
     102    c->cd(2);
     103    file.GetObject("MC/sparse/TrueEnergy/Impact", h2);
     104    h2->DrawCopy("colz");
     105    h2->DrawCopy("same");
     106    c->cd(3);
     107    file.GetObject("MC/theta/Impact", h2);
     108    h2->DrawCopy("colz");
     109    h2->DrawCopy("same");
    93110
    94111    // --------------------------------
  • trunk/FACT++/spectrum/simulation.sql

    r19915 r19971  
    1010      INTERVAL(LOG10(Energy), %108:sparse)  AS `.sparse_sim`,
    1111      INTERVAL(LOG10(Energy), %109:dense)  AS `.dense_sim`,
     12      INTERVAL(Impact/100, %110:impact)  AS `.impact`,
    1213
    13       (%110:spectrum)/POW(Energy, SpectralIndex) AS SpectralWeight,  -- FIXME: Is this correct for files with different Slopes?
     14      (%111:spectrum)/POW(Energy, SpectralIndex) AS SpectralWeight,  -- FIXME: Is this correct for files with different Slopes?
    1415      LogEnergyEst - log10(Energy) AS Residual
    1516   FROM
     
    2728   `.dense_est`,
    2829   `.dense_sim`,
     30   `.impact`,
    2931
    3032   -- Without any weight applied
     
    5961   ThetaDist USING(`.theta`)
    6062GROUP BY
    61    `.theta`, `.sparse_est`, `.sparse_sim`, `.dense_est`, `.dense_sim`
     63   `.theta`, `.sparse_est`, `.sparse_sim`, `.dense_est`, `.dense_sim`, `.impact`
  • trunk/FACT++/src/spectrum.cc

    r19913 r19971  
    139139    po::options_description binnings("Binnings");
    140140    binnings.add_options()
    141         ("theta",             var<Binning>(Binning(90, 0, 90)), "Add equidistant bins in theta (degrees). Syntax: N,lo,hi (Number N of equidistant bins between lo and hi)")
    142         ("theta-bin",         vars<double>(),                   "Add a bin-edge to the theta binning (degree)")
    143         ("energy-dense",      var<Binning>(Binning(30, 2, 5)),  "Add equidistant bins in log10 simulated energy. Syntax: N,lo,hi (Number N of equidistant bins between lo and hi)")
    144         ("energy-dense-bin",  vars<double>(),                   "Add a bin-edge to the binnnig in log10 simulated enegry")
    145         ("energy-sparse",     var<Binning>(Binning(15, 2, 5)),  "Add equidistant bins in log10 estimated energy. Syntax: N,lo,hi (Number N of equidistant bins between lo and hi)")
    146         ("energy-sparse-bin", vars<double>(),                   "Add a bin-edge to the binning in log10 estimated enegry")
     141        ("theta",             var<Binning>(Binning(90, 0, 90)),  "Add equidistant bins in theta (degrees). Syntax: N,lo,hi (Number N of equidistant bins between lo and hi)")
     142        ("theta-bin",         vars<double>(),                    "Add a bin-edge to the theta binning (degree)")
     143        ("energy-dense",      var<Binning>(Binning(30, 2, 5)),   "Add equidistant bins in log10 simulated energy. Syntax: N,lo,hi (Number N of equidistant bins between lo and hi)")
     144        ("energy-dense-bin",  vars<double>(),                    "Add a bin-edge to the binnnig in log10 simulated enegry")
     145        ("energy-sparse",     var<Binning>(Binning(15, 2, 5)),   "Add equidistant bins in log10 estimated energy. Syntax: N,lo,hi (Number N of equidistant bins between lo and hi)")
     146        ("energy-sparse-bin", vars<double>(),                    "Add a bin-edge to the binning in log10 estimated enegry")
     147        ("impact",            var<Binning>(Binning(28, 0, 280)), "Add equidistant bins in impact in meter. Syntax: N,lo,hi (Number N of equidistant bins between lo and hi)")
     148        ("impact-bin",        vars<double>(),                    "Add a bin-edge to the binning in impact in meter")
    147149        ;
    148150
     
    632634    Binning binning_dense  = conf.Get<Binning>("energy-dense");
    633635    Binning binning_sparse = conf.Get<Binning>("energy-sparse");
     636    Binning binning_impact = conf.Get<Binning>("impact");
    634637
    635638    cout << '\n';
    636     cout << "Binning 'theta':  " << binning_theta.str() << endl;
     639    cout << "Binning 'theta':  " << binning_theta.str()  << endl;
    637640    cout << "Binning 'dense':  " << binning_dense.str()  << endl;
    638     cout << "Binning 'sparse': " << binning_sparse.str()  << endl;
     641    cout << "Binning 'sparse': " << binning_sparse.str() << endl;
     642    cout << "Binning 'impact': " << binning_impact.str() << endl;
    639643
    640644    const uint16_t source_key = conf.Get<uint16_t>("source-key");
     
    656660    const string str_dense  = binning_dense.list();
    657661    const string str_sparse = binning_sparse.list();
     662    const string str_impact = binning_impact.list();
    658663
    659664    // -------------------------------------------------------------------------
     
    745750    CreateBinning(connection, qlog, binning_dense,  "Energy_dense",  "Dense binning in log10 Energy");
    746751    CreateBinning(connection, qlog, binning_sparse, "Energy_sparse", "Sparse binning in log10 Energy");
     752    CreateBinning(connection, qlog, binning_impact, "Impact",        "Binning in impact distance");
    747753
    748754    Dump(flog, connection, "BinningTheta");
    749755    Dump(flog, connection, "BinningEnergy_dense");
    750756    Dump(flog, connection, "BinningEnergy_sparse");
     757    Dump(flog, connection, "BinningImpact");
    751758
    752759    // -------------------------------------------------------------------
     
    928935        "   LEFT JOIN\n"
    929936        "      factmc.CorsikaSetup ON FileId=RUNNR\n"
    930         "   GROUP BY\n"
    931         "      `CSCAT[1]`, `CSCAT[2]`\n"
     937        // "   GROUP BY\n"
     938        // "      `CSCAT[1]`, `CSCAT[2]`\n"
    932939        "   ORDER BY\n"
    933940        "      MaxImpactHi\n"
     
    955962        cout << "Execution time: " << sec4 << "s\n\n";
    956963    }
    957 
    958964
    959965    // -------------------------------------------------------------------
     
    12971303        "   INDEX (`.dense_sim`)                 USING HASH,\n"
    12981304        "   INDEX (`.sparse_est`, `.sparse_sim`) USING HASH,\n"
    1299         "   INDEX (`.dense_est`, `.dense_sim`)   USING HASH\n"
     1305        "   INDEX (`.dense_est`, `.dense_sim`)   USING HASH,\n"
     1306        "   INDEX (`.impact`)                    USING HASH\n"
    13001307        ") ENGINE=MEMORY COMMENT='Sum of counts and (squared) weightes of Monte Carlo Data after analysis'\n"
    13011308        "AS\n"
     
    13171324
    13181325    //query6.template_defaults["columns"]   = "FileId, EvtNumber, CorsikaNumReuse,";
    1319     query8.template_defaults["columns"]   = "Energy, SpectralIndex,";
     1326    query8.template_defaults["columns"]   = "Energy, SpectralIndex, Impact,";
    13201327    query8.template_defaults["zenith"]    = "DEGREES(Theta)";
    13211328    query8.template_defaults["files"]     = "MonteCarloFiles";
     
    13271334    query8.template_defaults["dense"]     = str_dense.c_str();
    13281335    query8.template_defaults["theta"]     = str_theta.c_str();
     1336    query8.template_defaults["impact"]    = str_impact.c_str();
    13291337    query8.template_defaults["spectrum"]  = spectrum.c_str();
    13301338    query8.template_defaults["estimator"] = estimator.c_str();
     
    16191627        WriteHistogram(connection, hist_sim);
    16201628
     1629        // -------------------------------------------------------------------
     1630        // ------------------------ ImpactDistribution -----------------------
     1631        // -------------------------------------------------------------------
     1632
     1633        cout << separator("ImpactDistribution_"+*ib) << '\n';
     1634
     1635        Time start13;
     1636
     1637        mysqlpp::Query query13(&connection);
     1638        query13 <<
     1639            "CREATE TEMPORARY TABLE ImpactDistribution_%100:binning\n"
     1640            "(\n"
     1641            "   SignalN INT UNSIGNED NOT NULL\n"
     1642            ") ENGINE=MEMORY COMMENT='Impact Distribution of unweighted signal events after cuts'\n"
     1643            "AS\n"
     1644            "(\n"
     1645            "   SELECT\n"
     1646            "      %101:bin,\n"
     1647            "      `.impact`,\n"
     1648            "      SUM(SignalN) AS SignalN\n"
     1649            "   FROM\n"
     1650            "      AnalysisMC\n"
     1651            "   GROUP BY\n"
     1652            "      %101:bin, `.impact`\n"
     1653            "   ORDER BY\n"
     1654            "      %101:bin, `.impact`\n"
     1655            ")";
     1656
     1657        query13.parse();
     1658
     1659        query13.template_defaults["binning"] = ib->c_str();
     1660        query13.template_defaults["bin"]     = *ib=="theta" ? "`.theta`" : ("`."+*ib+"_sim`").c_str();
     1661
     1662        if (print_queries)
     1663            PrintQuery(query13.str());
     1664
     1665        qlog << query13 << ";\n" << endl;
     1666        if (connection.connected())
     1667        {
     1668            cout << query13.execute().info() << endl;
     1669            ShowWarnings(connection);
     1670            Dump(flog, connection, "ImpactDistribution_"+*ib);
     1671
     1672            const auto sec13 = Time().UnixTime()-start13.UnixTime();
     1673            cout << "Execution time: " << sec13 << "s\n";
     1674        }
     1675
     1676        hist_sim.name     = "Impact";
     1677        hist_sim.title    = "Distribution of Impact Distance";
     1678        hist_sim.y        = ".impact";
     1679        hist_sim.v        = "SignalN";
     1680        hist_sim.err      = "";
     1681        hist_sim.table    = "ImpactDistribution_"+*ib;
     1682        hist_sim.binningy = binning_impact;
     1683        hist_sim.axisy    = "Impact Distance [m]";
     1684        hist_sim.axisz    = "Counts";
     1685
     1686        WriteHistogram(connection, hist_sim);
     1687
     1688        // ===================================================================
     1689        // ===================================================================
     1690
    16211691        if (*ib=="theta")
    16221692            continue;
     1693
     1694        // ===================================================================
     1695        // ===================================================================
     1696
    16231697
    16241698        // -------------------------------------------------------------------
     
    17961870
    17971871        // -------------------------------------------------------------------
    1798         // ------------------------- SimulatedSpectrum -----------------------
     1872        // -------------------------- MigrationMatrix ------------------------
    17991873        // -------------------------------------------------------------------
    18001874
Note: See TracChangeset for help on using the changeset viewer.