Account for binwidth in HistogramModelFitter doctest.
authorW. Trevor King <wking@drexel.edu>
Fri, 5 Nov 2010 17:01:32 +0000 (13:01 -0400)
committerW. Trevor King <wking@drexel.edu>
Fri, 5 Nov 2010 17:01:32 +0000 (13:01 -0400)
Wider bins should increase the bin counts (all else being equal).

Also select a smaller sigma for the doctest, since that will make the
mean estimation less prone to errors large enough to break the test.

pysawsim/fit.py

index 2d35cf817a29ab95e2d470ba37d5ba1be126dce7..377d9835454fe05e28327d59c08eb8b902f6f482 100755 (executable)
@@ -315,7 +315,8 @@ class HistogramModelFitter (ModelFitter):
     ...         '''
     ...         p = params  # convenient alias
     ...         p[1] = abs(p[1])  # cannot have negative std. dev.
-    ...         self._model_data.counts = self.info['N']/(p[1]*sqrt(2*pi)) * (
+    ...         self._model_data.counts = (
+    ...             self.info['binwidth']*self.info['N']/(p[1]*sqrt(2*pi)) *
     ...             exp(-((self._model_data.bin_centers - p[0])/p[1])**2 / 2))
     ...         return self._model_data
     ...     def guess_initial_params(self):
@@ -330,8 +331,8 @@ class HistogramModelFitter (ModelFitter):
     ...             pstrings.append('%s=%.3f' % (name, param))
     ...         return ', '.join(pstrings)
 
-    >>> mu = 1.45
-    >>> sigma = 3.14
+    >>> mu = 3.14
+    >>> sigma = 1.45
     >>> data = array([gauss(mu,sigma) for i in range(int(1e4))])
     >>> h = Histogram()
     >>> h.from_data(data, h.calculate_bin_edges(data, 1))
@@ -340,7 +341,7 @@ class HistogramModelFitter (ModelFitter):
     >>> print m.model_string()
     p(x) = A exp((x-mu)^2 / (2 sigma^2))
     >>> print m.param_string(params)  # doctest: +ELLIPSIS
-    mu=1.4..., sigma=3.1...
+    mu=3.1..., sigma=1.4...
     """
     def set_data(self, data, info=None, rescale=False):
         self._data = data  # Histogram() instance
@@ -363,7 +364,7 @@ class HistogramModelFitter (ModelFitter):
         if self._rescale == True:
             params = [p*s for p,s in zip(params, self._param_scale_factors)]
         residual = self._data.counts - self.model(params).counts
-        if True:  # fit debugging
+        if False:  # fit debugging
             if not hasattr(self, '_i_'):
                 self._i_ = 0
             self._data.to_stream(open('hist-data.%d' % self._i_, 'w'))