Changeset 19904


Ignore:
Timestamp:
Dec 15, 2019, 4:23:56 PM (5 weeks ago)
Author:
tbretz
Message:
Implemented the Flux Sigma and control plots to check if a point should be replaced by an UL.
Location:
trunk/FACT++
Files:
3 edited

Legend:

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

    r19900 r19904  
    3232    fx.SetLineColor(kBlack);
    3333
    34     TF1 f0("f0", "0", 0, 90);
     34    TF1 f0("f0", "[0]", 0, 90);
     35    f0.SetParameter(0, 0);
    3536    f0.SetLineStyle(kDashed);
    3637    f0.SetLineColor(kBlack);
     
    211212
    212213    // --------------------------------
     214    f0.SetLineColor(kRed);
     215    // --------------------------------
    213216
    214217    c = new TCanvas("Integral Spectrum", "Integral Spectrum");
    215     c->SetLogy();
     218    c->Divide(2,2);
     219    c->cd(4);
     220    gPad->SetLogy();
    216221    file.GetObject("Data/Energy/Differential/IntegratedSpectrum", h1);
    217222    h1->SetLineColor(kGray);
    218223    h1->SetMarkerColor(kGray);
     224    if (h1->GetMinimum()<=0)
     225    {
     226        h1->SetMinimum(h1->GetMinimum(0)*0.1);
     227        cout << "WARNING: Integral Spectrum contains negative flux points!\n" << endl;
     228    }
    219229    h1->DrawCopy("P");
    220230    file.GetObject("Data/Energy/Integral/Spectrum", h1);
     
    238248    //h1->DrawCopy("P same");
    239249
     250    c->cd(1);
     251    gPad->SetGridx();
     252    file.GetObject("Data/Energy/Integral/SigmaFlux", h1);
     253    h1->SetMaximum(5);
     254    h1->SetMarkerStyle(kFullDotLarge);
     255    h1->DrawCopy("P");
     256    f0.SetParameter(0, 1);
     257    f0.DrawCopy("same");
     258
     259    c->cd(3);
     260    gPad->SetGridx();
     261    file.GetObject("Data/Energy/Integral/BackgroundI", h1);
     262    h1->Scale(5);
     263    h1->SetMaximum(30);
     264    h1->SetLineColor(kRed);
     265    h1->DrawCopy("P");
     266    file.GetObject("Data/Energy/Integral/SignalI", h1);
     267    h1->DrawCopy("P same");
     268    file.GetObject("Data/Energy/Integral/ExcessI", h1);
     269    h1->SetLineColor(kBlue);
     270    h1->DrawCopy("P same");
     271
     272    f0.SetParameter(0, 10);
     273    f0.DrawCopy("same");
     274
    240275    // --------------------------------
    241276
    242277    c = new TCanvas("Differential Spectrum", "Differential Spectrum");
    243     c->SetLogy();
     278    c->Divide(2,2);
     279    c->cd(4);
     280    gPad->SetLogy();
    244281    file.GetObject("Data/Energy/Differential/Spectrum", h1);
     282    if (h1->GetMinimum()<=0)
     283    {
     284        h1->SetMinimum(h1->GetMinimum(0)*0.1);
     285        cout << "WARNING: Differential Spectrum contains negative flux points!\n" << endl;
     286    }
    245287    h1->DrawCopy("P");
    246288
     
    261303    //h1->SetMarkerStyle(22);
    262304    //h1->DrawCopy("P same");
     305
     306    c->cd(1);
     307    gPad->SetGridx();
     308    file.GetObject("Data/Energy/Differential/SigmaFlux", h1);
     309    h1->SetMaximum(5);
     310    h1->SetMarkerStyle(kFullDotLarge);
     311    h1->DrawCopy("P");
     312    f0.SetParameter(0, 1);
     313    f0.DrawCopy("same");
     314
     315    c->cd(3);
     316    gPad->SetGridx();
     317    file.GetObject("Data/Energy/Differential/Background", h1);
     318    h1->Scale(5);
     319    h1->SetMaximum(30);
     320    h1->SetLineColor(kRed);
     321    h1->DrawCopy("P");
     322    file.GetObject("Data/Energy/Differential/Signal", h1);
     323    h1->DrawCopy("P same");
     324    file.GetObject("Data/Energy/Differential/Excess", h1);
     325    h1->SetLineColor(kBlue);
     326    h1->DrawCopy("P same");
     327
     328    f0.SetParameter(0, 10);
     329    f0.DrawCopy("same");
    263330}
  • trunk/FACT++/spectrum/spectrum.sql

    r19902 r19904  
    7474      %106:join2
    7575),
     76Flux AS
     77(
     78   SELECT -- Return final result
     79      *,
     80
     81      -- Differetial Spectrum
     82
     83      SimExcess/SimFluxW  AS  Efficiency,
     84
     85      Excess/SimExcess/Width/AreaTime  AS  ExcessRatio,
     86         1/SQRT(
     87             + POW(ErrExcess    / Excess,    2)
     88             + POW(ErrSimExcess / SimExcess, 2)
     89         )  AS  SigmaExcessRatio,
     90
     91
     92      Excess/SimExcess*SimFluxW/Width/AreaTime  AS  Flux,
     93         1/SQRT(
     94             + POW(ErrExcess    / Excess,    2)
     95             + POW(ErrSimExcess / SimExcess, 2)
     96             + POW(ErrSimFluxW  / SimFluxW,  2)
     97         )  AS  SigmaFlux,
     98
     99      -- Integral Spectrum
     100
     101      SimExcessI/SimFluxI  AS  EfficiencyI,
     102
     103      ExcessI/SimExcessI/AreaTime  AS  ExcessRatioI,
     104         1/SQRT(
     105             + POW(ErrExcessI    / ExcessI,    2)
     106             + POW(ErrSimExcessI / SimExcessI, 2)
     107         )  AS  SigmaExcessRatioI,
     108
     109
     110      ExcessI/SimExcessI*SimFluxI/AreaTime  AS  FluxI,
     111         1/SQRT(
     112             + POW(ErrExcessI    / ExcessI,    2)
     113             + POW(ErrSimExcessI / SimExcessI, 2)
     114             + POW(ErrSimFluxI   / SimFluxI,   2)
     115         )  AS  SigmaFluxI
     116
     117   FROM
     118      CombinedData
     119),
    76120Spectrum AS
    77121(
    78 SELECT -- Return final result
    79    *,
     122   SELECT -- Return final result
     123      *,
    80124
    81    -- Differetial Spectrum
    82 
    83    SimExcess/SimFluxW  AS  Efficiency,
    84 
    85    Excess/SimExcess/Width/AreaTime  AS  ExcessRatio,
    86    Excess/SimExcess/Width/AreaTime
    87       * SQRT(
    88           + POW(ErrExcess    / Excess,    2)
    89           + POW(ErrSimExcess / SimExcess, 2)
    90         )  AS  ErrExcessRatio,
    91 
    92 
    93    Excess/SimExcess*SimFluxW/Width/AreaTime  AS  Flux,
    94    Excess/SimExcess*SimFluxW/Width/AreaTime
    95       * SQRT(
    96           + POW(ErrExcess    / Excess,    2)
    97           + POW(ErrSimExcess / SimExcess, 2)
    98           + POW(ErrSimFluxW  / SimFluxW,  2)
    99         )  AS  ErrFlux,
    100 
    101    -- Integral Spectrum
    102 
    103    SimExcessI/SimFluxI  AS  EfficiencyI,
    104 
    105    ExcessI/SimExcessI/AreaTime  AS  ExcessRatioI,
    106    ExcessI/SimExcessI/AreaTime
    107       * SQRT(
    108           + POW(ErrExcessI    / ExcessI,    2)
    109           + POW(ErrSimExcessI / SimExcessI, 2)
    110         )  AS  ErrExcessRatioI,
    111 
    112 
    113    ExcessI/SimExcessI*SimFluxI/AreaTime  AS  FluxI,
    114    ExcessI/SimExcessI*SimFluxI/AreaTime
    115       * SQRT(
    116           + POW(ErrExcessI    / ExcessI,    2)
    117           + POW(ErrSimExcessI / SimExcessI, 2)
    118           + POW(ErrSimFluxI   / SimFluxI,   2)
    119         )  AS  ErrFluxI
    120 
    121 FROM
    122    CombinedData
    123 WHERE
    124    Excess>0
     125      ExcessRatio /SigmaExcessRatio   AS  ErrExcessRatio,
     126      ExcessRatioI/SigmaExcessRatioI  AS  ErrExcessRatioI,
     127      ABS(Flux)   /SigmaFlux          AS  ErrFlux,
     128      ABS(FluxI)  /SigmaFluxI         AS  ErrFluxI
     129   FROM
     130      Flux
    125131)
    126132SELECT
  • trunk/FACT++/src/spectrum.cc

    r19903 r19904  
    20602060            map<size_t, double> rolke_int;
    20612061
     2062            if (verbose>0)
     2063                cout << "Bin     Excess     Significance          Flux ErrFlux" << endl;
     2064
    20622065            for (auto ir=res13.cbegin(); ir!=res13.cend(); ir++)
    20632066            {
     
    20942097                if (verbose>0)
    20952098                {
    2096                     cout << setw(5)  << row["center"]  << ": ";
    2097                     cout << setw(10) << row["Excess"]  << " ";
    2098                     cout << setw(10) << row["Flux"]    << " ";
    2099                     cout << setw(10) << row["ErrFlux"] << " ";
    2100                     cout << setw(10) << row["Significance"] << '\n';
     2099                    cout << setw(5)  << row["center"] << ":";
     2100                    cout << " " << setw(10) << row["Excess"];
     2101                    cout << " " << setw(10) << row["Significance"];
     2102                    cout << " " << setw(10) << row["Flux"];
     2103                    cout << " " << setw(10) << row["ErrFlux"];
     2104                    cout << endl;
    21012105                }
    21022106            }
     
    21692173            WriteHistogram(connection, hist);
    21702174
     2175            hist.name = "SigmaFlux";
     2176            hist.v    = "SigmaFlux";
     2177            hist.err  = "";
     2178            WriteHistogram(connection, hist);
     2179
    21712180            if (*ib=="Energy")
    21722181            {
     
    21772186                hist.v     = "FluxI";
    21782187                hist.err   = "ErrFluxI";
     2188                WriteHistogram(connection, hist);
     2189
     2190                hist.dir   = "Data/Energy/Integral";
     2191                hist.name  = "SigmaFlux";
     2192                hist.v     = "SigmaFluxI";
     2193                hist.err   = "";
    21792194                WriteHistogram(connection, hist);
    21802195
Note: See TracChangeset for help on using the changeset viewer.