/* 

  Thermometer problem in Picat.

  From the Hakuru example thermometer.hk
  https://hakaru-dev.github.io/workflow/continuous/

  Cf my Gamble model gamble_thermometer.rkt

  This program was created by Hakan Kjellerstrand, hakank@gmail.com
  See also my Picat page: http://www.hakank.org/picat/

*/

import ppl_distributions, ppl_utils.
import util.

main => go.

/*

  var : m1
  mean = 20.9521
  HPD intervals:
  HPD interval (0.5): 17.00731139203729..24.96620524227995
  HPD interval (0.84): 12.09215541629798..29.36389337733698
  HPD interval (0.9): 10.36451003146743..30.82968066576282
  HPD interval (0.95): 8.45319336229872..33.23340342860692
  HPD interval (0.99): 3.79342020950252..37.47594851231587
  HPD interval (0.999): -1.63214465925364..42.36120922157517
  HPD interval (0.9999): -3.00393371835701..47.28067100088052
  HPD interval (0.99999): -3.00393371835701..51.13818397962348

  var : m2
  mean = 21.1086
  HPD intervals:
  HPD interval (0.5): 16.04259398150915..26.69301066375728
  HPD interval (0.84): 9.44234074716175..32.48807817435565
  HPD interval (0.9): 7.70961474132757..35.37547575694335
  HPD interval (0.95): 4.14576448712961..37.80303270040595
  HPD interval (0.99): -2.19554542341985..44.62767961781525
  HPD interval (0.999): -8.69771172667704..53.48819819580542
  HPD interval (0.9999): -16.11000597971537..56.58699050125793
  HPD interval (0.99999): -16.11000597971537..59.88694221132190

  var : noiseM
  mean = 2.49931
  HPD intervals:
  HPD interval (0.5): 1.63537940924772..3.10427227155504
  HPD interval (0.84): 1.38024641358305..3.87802211794910
  HPD interval (0.9): 1.01279150229543..3.69672951786627
  HPD interval (0.95): 1.05900856622449..3.89506068867401
  HPD interval (0.99): 1.01333480049546..3.98117851045969
  HPD interval (0.999): 1.00205708807430..3.99804672039954
  HPD interval (0.9999): 1.00034419866295..3.99877564143332
  HPD interval (0.99999): 1.00034419866295..3.99962386907992

  var : noiseT
  mean = 5.49944
  HPD intervals:
  HPD interval (0.5): 3.31906413394914..5.77056425938875
  HPD interval (0.84): 3.31906413394914..7.50170443602917
  HPD interval (0.9): 3.50366829638540..7.99497679294784
  HPD interval (0.95): 3.24275942484045..7.98673963127040
  HPD interval (0.99): 3.04174731906585..7.98732734470969
  HPD interval (0.999): 3.00595647376308..7.99960434855875
  HPD interval (0.9999): 3.00014402670792..7.99960434855875
  HPD interval (0.99999): 3.00014402670792..7.99989839969198

  var : t1
  mean = 20.9491
  HPD intervals:
  HPD interval (0.5): 17.41527400360535..24.50568985621916
  HPD interval (0.84): 13.13978079175492..28.66074725690494
  HPD interval (0.9): 11.92822317363790..30.40658780837546
  HPD interval (0.95): 9.48121361058946..32.23751452362503
  HPD interval (0.99): 5.34540398391366..36.25709226630411
  HPD interval (0.999): -0.13052910684177..39.79258917809992
  HPD interval (0.9999): -1.78996563662688..46.81711967752045
  HPD interval (0.99999): -1.78996563662688..49.01552880437639

  var : t2
  mean = 21.0821
  HPD intervals:
  HPD interval (0.5): 16.47568891029777..26.35849903841479
  HPD interval (0.84): 10.23849625593425..32.08940684965274
  HPD interval (0.9): 7.75224606251444..33.89807564899382
  HPD interval (0.95): 5.36399127610901..37.49798057711644
  HPD interval (0.99): -1.77533351859583..43.30425120955810
  HPD interval (0.999): -7.31205556151331..51.78777540734751
  HPD interval (0.9999): -10.72432788472920..59.21522033239848
  HPD interval (0.99999): -15.45655222294184..59.21522033239848



*/
go ?=>
  reset_store,
  run_model(10_000,$model,[mean,show_hpd_intervals]),
  nl,
  % fail,
  nl.
go => true.
  
  
model() =>
  NoiseT = uniform_dist(3,8),
  NoiseM = uniform_dist(1,4),

  T1 = normal_dist(21,NoiseT),
  T2 = normal_dist(T1,NoiseT),

  M1 = normal_dist(T1,NoiseM),
  M2 = normal_dist(T2,NoiseM),

  add("noiseT",NoiseT),
  add("noiseM",NoiseM),
  add("m1",M1),
  add("m2",M2),
  add("t1",T1),  
  add("t2",T2).