Thank you for the idea. I tried to use a Censored Poisson GLM with the following code:
with pm.Model() as poisson_glm_censored:
a = pm.HalfNormal('a', 10.0)
b = pm.Normal('b', 0.0, 1.0)
beta = lam * (0.5 * pm.math.erfc(a*x + b))
poisson_dist = pm.Poisson.dist(beta)
pm.Censored("y", poisson_dist, lower=0, upper=capacity, observed=y)
The prior simulations do look reasonable as far as I can see:

(blue dots = data points)
But the sampling completely fails, unfortunately:
with poisson_glm_censored:
trace_censored = pm.sample(tune=4000, step=[pm.NUTS(target_accept=0.9,max_treedepth=30)])
(I started with default parameters for pm.sample() and tried different sets of parameters all leading to bad sampling results.)
