already had a suspicion, that 66aea61 was wrong - did some tests and …

…rolled back to what is essentially the old algorithm + normalization for the time window input

playground/normalization_test.R

@@ -1,6 +1,26 @@
 library(ggplot2)
 library(magrittr)
 
+#### normalization experiments ####
+
+# normalization function
+nor <- function(x) { x/sum(x) }
+
+# test density vectors, all normalized
+a <- c(0.1,0.7,0.2)
+b <- c(0.3,0.4,0.3)
+c <- c(0.1,0.1,0.8)
+
+# sum -> all of these differ
+nor(a+b+c)
+nor(nor(a+b)+c)
+nor(a+nor(b+c))
+
+# product -> all are the same
+nor(a*b*c)
+nor(nor(a*b)*c)
+nor(a*nor(b*c))
+
 #### test simple sum ####
 
 system("currycarbon \"rangeBP(3000,2800)\" --densityFile /tmp/currycarbonOutput.tsv -q")
@@ -30,6 +50,19 @@ ggplot() +
     linewidth = 1, alpha = 0.5, color = "purple"
   )
 
+#### test two sums ####
+
+# this must not the same as the previous test
+system("currycarbon \"((3000,30) + (3000,30)) + rangeBP(3000,2800)\" --densityFile /tmp/currycarbonOutput.tsv -q")
+calPDFSum <- readr::read_tsv("/tmp/currycarbonOutput.tsv", col_types = readr::cols())
+
+ggplot() +
+  geom_line(
+    data = calPDFSum,
+    mapping = aes(x = yearBCAD, y = density), 
+    linewidth = 1, alpha = 0.5, color = "orange"
+  )
+
 #### test more complex arrangements ####
 
 system("currycarbon \"2900,100\" --densityFile /tmp/currycarbonOutput.tsv -q")
@@ -47,7 +80,7 @@ ggplot() +
     linewidth = 1, alpha = 0.5, color = "darkred"
   )
 
-system("currycarbon \"((3000,30) + rangeBP(3000,2800)) * (2900,100)\" --densityFile /tmp/currycarbonOutput.tsv -q")
+system("currycarbon \"(((3000,30) + (3000,30)) + rangeBP(3000,2800)) * (2900,100)\" --densityFile /tmp/currycarbonOutput.tsv -q")
 calPDFProd <- readr::read_tsv("/tmp/currycarbonOutput.tsv", col_types = readr::cols())
 
 ggplot() +
@@ -56,3 +89,5 @@ ggplot() +
     mapping = aes(x = yearBCAD, y = density), 
     linewidth = 1, alpha = 0.5, color = "darkgreen"
   )
+
+

src/Currycarbon/SumCalibration.hs

@@ -21,17 +21,21 @@ evalNamedCalExpr conf curve (NamedCalExpr exprID expr) =
 -- | Evaluate a dating expression by calibrating the individual dates and forming the respective
 --   sums and products of post-calibration density distributions
 evalCalExpr :: CalibrateDatesConf -> CalCurveBP -> CalExpr -> Either CurrycarbonException CalPDF
-evalCalExpr conf curve calExpr = normalize $ evalE calExpr
+evalCalExpr conf curve calExpr = norm $ evalE calExpr
     where
         evalE :: CalExpr -> Either CurrycarbonException CalPDF
-        evalE (UnCalDate a)    = normalize $ calibrateDate conf curve a
-        evalE (WindowBP a)     = normalize $ Right $ windowBP2CalPDF a
-        evalE (WindowBCAD a)   = normalize $ Right $ windowBCAD2CalPDF a
-        evalE (CalDate a)      = normalize $ Right a
-        evalE (SumCal a b)     = normalize $ eitherCombinePDFs (+) 0 (evalE a) (evalE b)
-        evalE (ProductCal a b) = normalize $ eitherCombinePDFs (*) 1 (evalE a) (evalE b)
-        normalize :: Either CurrycarbonException CalPDF -> Either CurrycarbonException CalPDF
-        normalize = mapEither id normalizeCalPDF
+        -- these are already normalized by their constructors
+        evalE (UnCalDate a)    = calibrateDate conf curve a
+        evalE (WindowBP a)     = Right $ windowBP2CalPDF a
+        evalE (WindowBCAD a)   = Right $ windowBCAD2CalPDF a
+        -- this can theoretically be non-normalized input
+        evalE (CalDate a)      = norm $ Right a
+        -- sums must not be normalized
+        evalE (SumCal a b)     = eitherCombinePDFs (+) 0 (evalE a) (evalE b)
+        -- products must be normalized (and their input, in case it's a sum)
+        evalE (ProductCal a b) = norm $ eitherCombinePDFs (*) 1 (norm $ evalE a) (norm $ evalE b)
+        norm :: Either CurrycarbonException CalPDF -> Either CurrycarbonException CalPDF
+        norm = mapEither id normalizeCalPDF
         -- https://hackage.haskell.org/package/either-5.0.2/docs/Data-Either-Combinators.html
         mapEither :: (a -> c) -> (b -> d) -> Either a b -> Either c d
         mapEither f _ (Left x)  = Left (f x)
@@ -84,7 +88,7 @@ windowBCAD2CalPDF :: TimeWindowBCAD -> CalPDF
 windowBCAD2CalPDF (TimeWindowBCAD name start stop) =
     let years = VU.fromList $ [start..stop]
         dens = VU.replicate (VU.length years) 1
-    in CalPDF name years dens
+    in normalizeCalPDF $ CalPDF name years dens
 
 windowBP2CalPDF :: TimeWindowBP -> CalPDF
 windowBP2CalPDF (TimeWindowBP name start stop) =