summaryrefslogtreecommitdiff
path: root/farming/statistics.lua
diff options
context:
space:
mode:
authorVanessa Ezekowitz <vanessaezekowitz@gmail.com>2017-02-19 12:01:09 -0500
committerVanessa Ezekowitz <vanessaezekowitz@gmail.com>2017-02-19 12:01:09 -0500
commita7bebc671920fd58992de4ff5a5f2716c891565c (patch)
tree33fdcbfefeffe292e1c20c156d84cfafb89b8b80 /farming/statistics.lua
parent2aae75d105c39af6fa2b5a84672726ee5c19fec1 (diff)
downloaddreambuilder_modpack-a7bebc671920fd58992de4ff5a5f2716c891565c.tar
dreambuilder_modpack-a7bebc671920fd58992de4ff5a5f2716c891565c.tar.gz
dreambuilder_modpack-a7bebc671920fd58992de4ff5a5f2716c891565c.tar.bz2
dreambuilder_modpack-a7bebc671920fd58992de4ff5a5f2716c891565c.tar.xz
dreambuilder_modpack-a7bebc671920fd58992de4ff5a5f2716c891565c.zip
updated bobblocks, gloopblocks, maptools, mesecons, moreblocks,
moreores, plantlife, and xban2, and switched to TenPlus1's "farming redo" mod in lieu of Pilzadam's farming_plus since the latter breaks default farming.
Diffstat (limited to 'farming/statistics.lua')
-rw-r--r--farming/statistics.lua174
1 files changed, 174 insertions, 0 deletions
diff --git a/farming/statistics.lua b/farming/statistics.lua
new file mode 100644
index 0000000..c8defa2
--- /dev/null
+++ b/farming/statistics.lua
@@ -0,0 +1,174 @@
+local statistics = {}
+local ROOT_2 = math.sqrt(2.0)
+
+-- Approximations for erf(x) and erfInv(x) from
+-- https://en.wikipedia.org/wiki/Error_function
+
+local erf
+local erf_inv
+
+local A = 8 * (math.pi - 3.0)/(3.0 * math.pi * (4.0 - math.pi))
+local B = 4.0 / math.pi
+local C = 2.0/(math.pi * A)
+local D = 1.0 / A
+
+erf = function(x)
+
+ if x == 0 then return 0; end
+
+ local xSq = x * x
+ local aXSq = A * xSq
+ local v = math.sqrt(1.0 - math.exp(-xSq * (B + aXSq) / (1.0 + aXSq)))
+
+ return (x > 0 and v) or -v
+end
+
+erf_inv = function(x)
+
+ if x == 0 then return 0; end
+
+ if x <= -1 or x >= 1 then return nil; end
+
+ local y = math.log(1 - x * x)
+ local u = C + 0.5 * y
+ local v = math.sqrt(math.sqrt(u * u - D * y) - u)
+
+ return (x > 0 and v) or -v
+end
+
+local function std_normal(u)
+ return ROOT_2 * erf_inv(2.0 * u - 1.0)
+end
+
+local poisson
+local cdf_table = {}
+
+local function generate_cdf(lambda_index, lambda)
+
+ local max = math.ceil(4 * lambda)
+ local pdf = math.exp(-lambda)
+ local cdf = pdf
+ local t = { [0] = pdf }
+
+ for i = 1, max - 1 do
+ pdf = pdf * lambda / i
+ cdf = cdf + pdf
+ t[i] = cdf
+ end
+
+ return t
+end
+
+for li = 1, 100 do
+ cdf_table[li] = generate_cdf(li, 0.25 * li)
+end
+
+poisson = function(lambda, max)
+
+ if max < 2 then
+ return (math.random() < math.exp(-lambda) and 0) or 1
+ elseif lambda >= 2 * max then
+ return max
+ end
+
+ local u = math.random()
+ local lambda_index = math.floor(4 * lambda + 0.5)
+ local cdfs = cdf_table[lambda_index]
+
+ if cdfs then
+
+ lambda = 0.25 * lambda_index
+
+ if u < cdfs[0] then return 0; end
+ if max > #cdfs then max = #cdfs + 1 else max = math.floor(max); end
+ if u >= cdfs[max - 1] then return max; end
+
+ if max > 4 then -- Binary search
+
+ local s = 0
+
+ while s + 1 < max do
+
+ local m = math.floor(0.5 * (s + max))
+
+ if u < cdfs[m] then max = m; else s = m; end
+ end
+ else
+ for i = 1, max - 1 do
+ if u < cdfs[i] then return i; end
+ end
+ end
+
+ return max
+ else
+ local x = lambda + math.sqrt(lambda) * std_normal(u)
+
+ return (x < 0.5 and 0) or (x >= max - 0.5 and max) or math.floor(x + 0.5)
+ end
+end
+
+-- Error function.
+statistics.erf = erf
+
+-- Inverse error function.
+statistics.erf_inv = erf_inv
+
+--- Standard normal distribution function (mean 0, standard deviation 1).
+ --
+ -- @return
+ -- Any real number (actually between -3.0 and 3.0).
+
+statistics.std_normal = function()
+
+ local u = math.random()
+
+ if u < 0.001 then
+ return -3.0
+ elseif u > 0.999 then
+ return 3.0
+ end
+
+ return std_normal(u)
+end
+
+--- Standard normal distribution function (mean 0, standard deviation 1).
+ --
+ -- @param mu
+ -- The distribution mean.
+ -- @param sigma
+ -- The distribution standard deviation.
+ -- @return
+ -- Any real number (actually between -3*sigma and 3*sigma).
+
+statistics.normal = function(mu, sigma)
+
+ local u = math.random()
+
+ if u < 0.001 then
+ return mu - 3.0 * sigma
+ elseif u > 0.999 then
+ return mu + 3.0 * sigma
+ end
+
+ return mu + sigma * std_normal(u)
+end
+
+--- Poisson distribution function.
+ --
+ -- @param lambda
+ -- The distribution mean and variance.
+ -- @param max
+ -- The distribution maximum.
+ -- @return
+ -- An integer between 0 and max (both inclusive).
+
+statistics.poisson = function(lambda, max)
+
+ lambda, max = tonumber(lambda), tonumber(max)
+
+ if not lambda or not max or lambda <= 0 or max < 1 then return 0; end
+
+ return poisson(lambda, max)
+end
+
+return statistics