code: def _call(self, x):
          
          y = self.data
          r = self.background
  
          obj = self.domain.zero()
  
          # x + r - y + y * log(y / (x + r)) = x - y * log(x + r) + c1
          # with c1 = r - y + y * log y
          i = x.ufuncs.greater_equal(0)
          obj[i] = x[i] + r[i] - y[i]
  
          j = y.ufuncs.greater(0)
          k = i.ufuncs.logical_and(j)
          obj[k] += y[k] * (y[k] / (x[k] + r[k])).ufuncs.log()
  
          # (y / (2 * r^2)) * x^2 + (1 - y / r) * x + r - b + b * log(b / r)
          # = (y / (2 * r^2)) * x^2 + (1 - y / r) * x + c2
          # with c2 = r - b + b * log(b / r)
          i = i.ufuncs.logical_not()
          obj[i] += (y[i] / (2 * r[i]**2) * x[i]**2 + (1 - y[i] / r[i]) * x[i] +
                     r[i] - y[i])
  
          k = i.ufuncs.logical_and(j)
          obj[k] += y[k] * (y[k] / r[k]).ufuncs.log()
  
          return obj.inner(self.domain.one())
docs: Return the KL-diveregnce in the point ``x``.
  
          If any components of ``x`` is non-positive, the value is positive
          infinity.