Index: trunk/FACT++/spectrum/analysis.sql =================================================================== --- trunk/FACT++/spectrum/analysis.sql (revision 19890) +++ trunk/FACT++/spectrum/analysis.sql (revision 19893) @@ -8,11 +8,13 @@ - `LogEnergyEst` logarithm of estimated energy in GeV - In additon, all columns provided by the %0:column clause must be + In additon, all columns provided by the 100-clause must be returned. (Note that you must not add a comma behind it) - %101:files table containing the `FileId`s to analyze. - %102:runinfo table with the run info data - %103:events table with the image parameters - %104:positions table with the source positions in the camera + 100| %100:columns + 101| %101:files:: table containing the `FileId`s to analyze. + 102| %102:runinfo:: table with the run info data + 103| %103:events:: table with the image parameters + 104| %104:positions:: table with the source positions in the camera + 105| %105:estimator:: estimator for log10 energy *************************************************************************** */ Index: trunk/FACT++/spectrum/data.sql =================================================================== --- trunk/FACT++/spectrum/data.sql (revision 19890) +++ trunk/FACT++/spectrum/data.sql (revision 19893) @@ -1,25 +1,17 @@ -WITH Table0 AS -( - SELECT - INTERVAL(LogEnergyEst, %106:bins) AS `.energy`, - COUNT(IF(Weight>0, 1, NULL)) AS `Signal`, - COUNT(IF(Weight<0, 1, NULL))/5 AS `Background` - -- SUM(IF(Weight>0, 1./log10(Size), 0)) AS `Signal`, - -- SUM(IF(Weight<0, 1./log10(Size), 0)) AS `Background` +SELECT + INTERVAL(fZenithDistanceMean, %106:theta) AS `.theta`, + INTERVAL(LogEnergyEst, %107:sparse) AS `.sparse_est`, - -- FIXME: Add excess vs theta - FROM - Excess - GROUP BY - `.energy` - ORDER BY - `.energy` -) + -- Signal and Background counts + COUNT(IF(Weight>0, 1, NULL)) AS `Signal`, + COUNT(IF(Weight<0, 1, NULL)) AS `Background`, -SELECT - *, - `Signal` - `Background` AS `Excess`, - LiMa(`Signal`, `Background`) AS `Significance`, - ExcErr(`Signal`, `Background`) AS `ErrExcess` + -- Average Energy + SUM(Weight*POW(10, LogEnergyEst)) AS `SumEnergyEst`, + SUM(Weight) AS `SumW` FROM - Table0 + Excess +GROUP BY + `.theta`, `.sparse_est` +ORDER BY + `.theta`, `.sparse_est` Index: trunk/FACT++/spectrum/simulation.sql =================================================================== --- trunk/FACT++/spectrum/simulation.sql (revision 19890) +++ trunk/FACT++/spectrum/simulation.sql (revision 19893) @@ -2,10 +2,12 @@ ( SELECT - Weight, - INTERVAL(Zd, %106:theta) AS `.theta`, - INTERVAL(LogEnergyEst, %107:energyest) AS `.energyest`, - INTERVAL(log10(Energy), %108:energysim) AS `.energysim`, - (%109:spectrum)/pow(Energy, SpectralIndex) AS SpectralWeight, -- FIXME: Is this correct for files with different Slopes? + Weight, Energy, LogEnergyEst, + INTERVAL(Zd, %106:theta) AS `.theta`, + INTERVAL(LogEnergyEst, %107:sparse) AS `.sparse_est`, + INTERVAL(LogEnergyEst, %108:dense) AS `.dense_est`, + INTERVAL(LOG10(Energy), %107:sparse) AS `.sparse_sim`, + INTERVAL(LOG10(Energy), %108:dense) AS `.dense_sim`, + (%109:spectrum)/POW(Energy, SpectralIndex) AS SpectralWeight, -- FIXME: Is this correct for files with different Slopes? LogEnergyEst - log10(Energy) AS Residual FROM @@ -16,72 +18,41 @@ -- INNER JOIN -- factmc.RunInfoMC USING(FIleId) -), +) +SELECT + `.theta`, + `.sparse_est`, + `.sparse_sim`, + `.dense_est`, + `.dense_sim`, -Table1 AS -( - SELECT - `.energyest`, - `.energysim`, + -- Without any weight applied + COUNT(IF(Weight>0, 1, NULL)) AS SignalN, + COUNT(IF(Weight<0, 1, NULL)) AS BackgroundN, - -- Signal, Background, Excess - SUM( IF(Weight>0, ZdWeight*SpectralWeight, 0)) OVER EnergyEst AS `SignalW`, - SUM( IF(Weight<0, ZdWeight*SpectralWeight, 0)) OVER EnergyEst AS `BackgroundW`, - SUM( IF(Weight>0, POW(ErrZdWeight*SpectralWeight,2), 0)) OVER EnergyEst AS `SignalW2`, - SUM( IF(Weight<0, POW(ErrZdWeight*SpectralWeight,2), 0)) OVER EnergyEst AS `BackgroundW2`, - COUNT(IF(Weight>0, 1, NULL)) OVER EnergyEst AS `SignalN`, - COUNT(IF(Weight<0, 1, NULL)) OVER EnergyEst AS `BackgroundN`, - -- SUM( IF(Table9.Weight>0, ThetaHist.Weight*spectrum/Energy, 0)) OVER EnergyEst AS `SignalW`, - -- SUM( IF(Table9.Weight<0, ThetaHist.Weight*spectrum/Energy, 0)) OVER EnergyEst AS `BackgroundW`, - -- SUM( IF(Table9.Weight>0, POW(ThetaHist.Weight*spectrum/Energy,2), 0)) OVER EnergyEst AS `SignalW2`, - -- SUM( IF(Table9.Weight<0, POW(ThetaHist.Weight*spectrum/Energy,2), 0)) OVER EnergyEst AS `BackgroundW2`, + -- Without ZdWeight applied +/* + SUM( IF(Weight>0, SpectralWeight, 0)) AS Signal, + SUM( IF(Weight<0, SpectralWeight, 0)) AS Background, + SUM( IF(Weight>0, POW(SpectralWeight,2), 0)) AS Signal2, + SUM( IF(Weight<0, POW(SpectralWeight,2), 0)) AS Background2, +*/ - -- Threshold - SUM( Weight * ZdWeight*SpectralWeight ) OVER EnergySim AS `ThresholdW`, - SUM(POW(Weight * ErrZdWeight*SpectralWeight,2)) OVER EnergySim AS `ThresholdW2`, - SUM( Weight ) OVER EnergySim AS `ThresholdN`, - -- SUM( Table9.Weight/Energy * ThetaHist.Weight*spectrum ) OVER EnergySim AS `ThresholdW`, - -- SUM( POW(Table9.Weight/Energy * ThetaHist.Weight*spectrum,2) ) OVER EnergySim AS `ThresholdW2`, - -- SUM( Table9.Weight/Energy ) OVER EnergySim AS `ThresholdN` + -- Binning in estimated energy: Signal, Background, Excess + SUM( IF(Weight>0, ZdWeight*SpectralWeight, 0)) AS SignalW, + SUM( IF(Weight<0, ZdWeight*SpectralWeight, 0)) AS BackgroundW, + SUM( IF(Weight>0, POW(ErrZdWeight*SpectralWeight,2), 0)) AS SignalW2, + SUM( IF(Weight<0, POW(ErrZdWeight*SpectralWeight,2), 0)) AS BackgroundW2, - -- Estimators - SUM(IF(Weight>0, ZdWeight*SpectralWeight, 0)) OVER EnergySim AS SimW, - SUM(IF(Weight>0, Residual* ZdWeight*SpectralWeight, 0)) OVER EnergyEst AS EstSum, - SUM(IF(Weight>0, Residual* ZdWeight*SpectralWeight, 0)) OVER EnergySim AS SimSum, - SUM(IF(Weight>0, POW(Residual,2)*ZdWeight*SpectralWeight, 0)) OVER EnergyEst AS EstSum2, - SUM(IF(Weight>0, POW(Residual,2)*ZdWeight*SpectralWeight, 0)) OVER EnergySim AS SimSum2, + -- Energy Estimation + SUM(IF(Weight>0, Residual* ZdWeight*SpectralWeight, 0)) AS ResidualW, + SUM(IF(Weight>0, POW(Residual,2)*ZdWeight*SpectralWeight, 0)) AS ResidualW2, - -- Migration - SUM(Weight * ZdWeight*SpectralWeight) OVER Migration AS `MigrationW`, - SUM(Weight ) OVER Migration AS `MigrationN` - - -- FIXME: Add ExcessN vs theta - -- FIXME: Add ExcessW vs theta - - FROM - Table0 - INNER JOIN - ThetaDistribution USING(`.theta`) - WINDOW - EnergyEst AS (PARTITION BY `.energyest`), - EnergySim AS (PARTITION BY `.energysim`), - Migration AS (PARTITION BY `.energysim`,`.energyest`) -) - -SELECT DISTINCT - - *, - - `SignalW` - `BackgroundW`/5 AS `ExcessW`, - `SignalN` - `BackgroundN`/5 AS `ExcessN`, - - ExcErr(`SignalW2`, `BackgroundW2`/5) AS `ErrExcessW`, - ExcErr(`SignalN`, `BackgroundN` /5) AS `ErrExcessN`, - - IF(SignalW=0, 0, EstSum / SignalW) AS BiasEst, -- FIMXE: Is NULL better? - IF(SimW =0, 0, SimSum / SimW) AS BiasSim, -- FIMXE: Is NULL better? - - IF(SignalW=0, 0, SQRT(EstSum2/SignalW - POW(EstSum/SignalW, 2))) AS ResolutionEst, - IF(SimW =0, 0, SQRT(SimSum2/SimW - POW(SimSum/SimW, 2))) AS ResolutionSim + SUM(IF(Weight>0, Energy *ZdWeight*SpectralWeight, 0)) AS SumEnergySimW, + SUM(IF(Weight>0, POW(10, LogEnergyEst)*ZdWeight*SpectralWeight, 0)) AS SumEnergyEstW FROM - Table1 + Table0 +INNER JOIN + ThetaDist USING(`.theta`) +GROUP BY + `.theta`, `.sparse_est`, `.sparse_sim`, `.dense_est`, `.dense_sim` Index: trunk/FACT++/spectrum/spectrum.sql =================================================================== --- trunk/FACT++/spectrum/spectrum.sql (revision 19893) +++ trunk/FACT++/spectrum/spectrum.sql (revision 19893) @@ -0,0 +1,74 @@ +WITH BinnedData AS -- Bin data and create sums per bin +( + SELECT + %101:bin, + SUM(`Signal`) AS `Signal`, + SUM(`Background`)/5 AS `Background`, + SUM(SumEnergyEst) AS SumEnergyEst, + SUM(SumW) AS SumW + FROM + SummaryData + GROUP BY + %101:bin +), +ResultData AS -- Calculate results for all bins +( + SELECT + *, + SQRT(`Signal`) AS ErrSignal, + SQRT(`Background`) AS ErrBackground, + `Signal` - `Background` AS Excess, + ExcErr(`Signal`, `Background`) AS ErrExcess, + LiMa( `Signal`, `Background`) AS Significance + FROM + BinnedData +), +CombinedData AS -- Joind data together and for conveninece (easier reading of the query) rename the columns +( + SELECT + %102:id.%101:bin, + %102:id.lo, + %102:id.hi, + (%102:id.lo+%102:id.hi)/2 AS center, + %102:id.Scale, + Data.Excess, + Data.Signal, + Data.Background, + Data.ErrExcess, + Data.ErrSignal, + Data.ErrBackground, + Significance, + SumEnergyEst/SumW AS AvgEnergyEst, + %102:id.ExcessFluxW AS SimExcess, + %102:id.ErrExcessFluxW AS ErrSimExcess, + Sim.SimFluxW/%103:weight AS SimFluxW, -- Possible correcton for already applied ZdWeights + Sim.ErrSimFluxW/%104:errweight AS ErrSimFluxW -- Possible correcton for already applied ZdWeights + FROM + ResultData Data + INNER JOIN + %105:join1 + INNER JOIN + %106:join2 +) +SELECT -- Return final result + *, + SimExcess/SimFluxW AS Efficiency, + Excess/SimExcess*Scale AS ExcessRatio, + Excess/SimExcess*Scale + * SQRT( + + POW(ErrExcess / Excess, 2) + + POW(ErrSimExcess / SimExcess, 2) + ) AS ErrExcessRatio, + Excess/SimExcess*SimFluxW*Scale AS Flux, + Excess/SimExcess*SimFluxW*Scale + * SQRT( + + POW(ErrExcess / Excess, 2) + + POW(ErrSimExcess / SimExcess, 2) + + POW(ErrSimFluxW / SimFluxW, 2) + ) AS ErrFlux +FROM + CombinedData +WHERE + Excess>0 +ORDER BY + %101:bin Index: trunk/FACT++/spectrum/summary-est.sql =================================================================== --- trunk/FACT++/spectrum/summary-est.sql (revision 19893) +++ trunk/FACT++/spectrum/summary-est.sql (revision 19893) @@ -0,0 +1,58 @@ +WITH Theta AS -- Get Total Observation Time +( + SELECT SUM(OnTime) AS ObsTime FROM ThetaDist +), +Area AS -- Get total simulated area +( + SELECT POW(MinImpactHi,2)*PI() AS Area FROM MonteCarloArea +), +Binning AS -- Get binning and calculate weights for each bin +( + SELECT + *, + 1000/(POW(10,hi)-POW(10,lo))/Area/ObsTime AS Scale + FROM + BinningEnergy_%100:binning + CROSS JOIN + Theta, Area +), +Analyzed AS -- Summarizy data after cuts (analyzed) in bins +( + SELECT + `.%100:binning:_est`, + SUM(SignalN) AS SignalN, + SUM(SignalW) AS SignalW, + SUM(SignalW2) AS SignalW2, + SUM(BackgroundN) AS BackgroundN, + SUM(BackgroundW) AS BackgroundW, + SUM(BackgroundW2) AS BackgroundW2, + SUM(ResidualW) AS ResidualW, + SUM(ResidualW2) AS ResidualW2, + SUM(SumEnergySimW) AS SumEnergySimW, + SUM(SumEnergyEstW) AS SumEnergyEstW + FROM + AnalysisMC + GROUP BY + `.%100:binning:_est` +) +SELECT -- Return Result + *, + SignalN - BackgroundN/5 AS ExcessN, + SQRT(SignalN + BackgroundN/25) AS ErrExcessN, + SQRT(SignalN) AS ErrSignalN, + SQRT(BackgroundN) AS ErrBackgroundN, + Scale*SignalW AS SignalFluxW, + Scale*BackgroundW/5 AS BackgroundFluxW, + Scale*(SignalW - BackgroundW/5) AS ExcessFluxW, + Scale*SQRT(SignalW2) AS ErrSignalFluxW, + Scale*SQRT(BackgroundW2)/5 AS ErrBackgroundFluxW, + Scale*SQRT(SignalW2 + BackgroundW2/25) AS ErrExcessFluxW, + SumEnergyEstW/SignalW AS AvgEnergyEstW, + SumEnergySimW/SignalW AS AvgEnergySimW, + IF(SignalW>0, ResidualW/SignalW, NULL) AS BiasW, + IF(SignalW>0, ResidualW/SignalW*SQRT(ResidualW2/POW(ResidualW,2) + SignalW2/POW(SignalW,2)), NULL) AS ErrBiasW, + IF(SignalW>0, SQRT(ResidualW2/SignalW - POW(ResidualW/SignalW, 2)), NULL) AS ResolutionW +FROM + Analyzed +INNER JOIN + Binning ON `.%100:binning:_est`=bin Index: trunk/FACT++/spectrum/summary-sim.sql =================================================================== --- trunk/FACT++/spectrum/summary-sim.sql (revision 19893) +++ trunk/FACT++/spectrum/summary-sim.sql (revision 19893) @@ -0,0 +1,117 @@ +WITH Theta AS -- Get total Observation time +( + SELECT SUM(OnTime) AS ObsTime FROM ThetaDist +), +Area AS -- Get total simulated area +( + SELECT POW(MinImpactHi,2)*PI() AS Area FROM MonteCarloArea +), +Binning AS -- Create binning and calculate weights per bin +( + SELECT + *, + (%103:binwidth)/Area/ObsTime AS Scale + FROM + %101:binning + CROSS JOIN + Theta, Area +), +Simulated AS -- Summarize simulated (corsika) data in bins +( + SELECT + %102:bin, + SUM(CountN) AS SimCountN, + SUM(SumW) AS SimSumW, + SUM(SumW2) AS SimSumW2 + FROM + WeightedOriginalMC + GROUP BY + %102:bin +), +Triggered AS -- Summarizy triggered (ceres) data in bins +( + SELECT + %102:bin, + SUM(CountN) AS TrigCountN, + SUM(SumW) AS TrigSumW, + SUM(SumW2) AS TrigSumW2 + FROM + WeightedEventsMC + GROUP BY + %102:bin +), +Analyzed AS -- Summarize data after cuts (analyzed) in bins +( + SELECT + %102:bin, + SUM(SignalN) AS SignalN, + SUM(SignalW) AS SignalW, + SUM(SignalW2) AS SignalW2, + SUM(BackgroundN) AS BackgroundN, + SUM(BackgroundW) AS BackgroundW, + SUM(BackgroundW2) AS BackgroundW2, + SUM(ResidualW) AS ResidualW, + SUM(ResidualW2) AS ResidualW2, + SUM(SumEnergySimW) AS SumEnergyEstW, + SUM(SumEnergyEstW) AS SumEnergySimW + FROM + AnalysisMC + GROUP BY + %102:bin +), +CombinedData AS -- Combine all Data together +( + SELECT + *, + TrigCountN/SimCountN AS TriggerEfficiencyN, + TrigSumW/SimSumW AS TriggerEfficiencyW, + SignalW - BackgroundW/5 AS ExcessW, + SignalW2 + BackgroundW2/25 AS ExcessW2 + FROM + Simulated + INNER JOIN + Triggered USING (%102:bin) + INNER JOIN + Analyzed USING (%102:bin) + INNER JOIN + Binning ON Simulated.%102:bin=bin +), +Table0 AS -- Derive valued +( + SELECT -- Everything scaled is "Flux", Everything unscaled is without "Flux", Corsika Data is Sim*, Triggered Data is Trig* and data after cuts is without prefix + *, + SQRT(SimCountN) AS ErrSimCountN, + SQRT(TrigCountN) AS ErrTrigCountN, + SQRT(SignalN) AS ErrSignalN, + Scale*SimSumW AS SimFluxW, + Scale*TrigSumW AS TrigFluxW, + Scale*SignalW AS SignalFluxW, + Scale*BackgroundW/5 AS BackgroundFluxW, + Scale*ExcessW AS ExcessFluxW, + Scale*SQRT(SimSumW2) AS ErrSimFluxW, + Scale*SQRT(TrigSumW2) AS ErrTrigFluxW, + Scale*SQRT(SignalW2) AS ErrSignalFluxW, + Scale*SQRT(BackgroundW2)/5 AS ErrBackgroundFluxW, + Scale*SQRT(ExcessW2) AS ErrExcessFluxW, + SumEnergyEstW/SignalW AS AvgEnergyEstW, + SumEnergySimW/SignalW AS AvgEnergySimW, + Area*TriggerEfficiencyN AS EffectiveAreaN, + Area*TriggerEfficiencyW AS EffectiveAreaW, + SignalN/TrigCountN AS CutEfficiencyN, + IF(ExcessW<0, 0, ExcessW/TrigSumW) AS CutEfficiencyW, + TriggerEfficiencyN * SQRT(1/TrigCountN + 1/SimCountN) AS ErrTriggerEfficiencyN, + TriggerEfficiencyW * SQRT(TrigSumW2/POW(TrigSumW,2) + SimSumW2/POW(SimSumW,2)) AS ErrTriggerEfficiencyW, + IF(SignalW>0, ResidualW/SignalW, NULL) AS BiasW, + IF(SignalW>0, ResidualW/SignalW*SQRT(ResidualW2/POW(ResidualW,2) + SignalW2/POW(SignalW,2)), NULL) AS ErrBiasW, + IF(SignalW>0, SQRT(ResidualW2/SignalW - POW(ResidualW/SignalW, 2)), NULL) AS ResolutionW + FROM + CombinedData +) +SELECT -- Return derived values and result + *, + Area*ErrTriggerEfficiencyN AS ErrEffectiveAreaN, + Area*ErrTriggerEfficiencyW AS ErrEffectiveAreaW, + CutEfficiencyN*SQRT(1/SignalN + 1/TrigCountN) AS ErrCutEfficiencyN, + CutEfficiencyW*SQRT(ExcessW2/POW(ExcessW,2) + TrigSumW2/POW(TrigSumW,2)) AS ErrCutEfficiencyW +FROM + Table0