docs: def pdm (t, x, u, periods, nbin, nshift=8, nsmc=256, numc=256, weights=False, parallel=True):
      """Perform phase dispersion minimization.
  
      t : 1D array
        time coordinate
      x : 1D array, same size as *t*
        observed value
      u : 1D array, same size as *t*
        uncertainty on observed value; same units as `x`
      periods : 1D array
        set of candidate periods to sample; same units as `t`
      nbin : int
        number of phase bins to construct
      nshift : int=8
        number of shifted binnings to sample to combact statistical flukes
      nsmc : int=256
        number of Monte Carlo shufflings to compute, to evaluate the
        significance of the minimal theta value.
      numc : int=256
        number of Monte Carlo added-noise datasets to compute, to evaluate
        the uncertainty in the location of the minimal theta value.
      weights : bool=False
        if True, 'u' is actually weights, not uncertainties.
        Usually weights = u**-2.
      parallel : default True
        Controls parallelization of the algorithm. Default
        uses all available cores. See `pwkit.parallel.make_parallel_helper`.
  
      Returns named tuple of:
  
      thetas : 1D array
        values of theta statistic, same size as `periods`
      imin
        index of smallest (best) value in `thetas`
      pmin
        the `period` value with the smallest (best) `theta`
      mc_tmins
        1D array of size `nsmc` with Monte Carlo samplings of minimal
        theta values for shufflings of the data; assesses significance of the peak
      mc_pvalue
        probability (between 0 and 1) of obtaining the best theta value
        in a randomly-shuffled dataset
      mc_pmins
        1D array of size `numc` with Monte Carlo samplings of best
        period values for noise-added data; assesses uncertainty of `pmin`
      mc_puncert
        standard deviation of `mc_pmins`; approximate uncertainty
        on `pmin`.
  
      We don't do anything clever, so runtime scales at least as
      ``t.size * periods.size * nbin * nshift * (nsmc + numc + 1)``.
  
      """
      t = np.asfarray (t)
      x = np.asfarray (x)
      u = np.asfarray (u)
      periods = np.asfarray (periods)
      t, x, u, periods = np.atleast_1d (t, x, u, periods)
      nbin = int (nbin)
      nshift = int (nshift)
      nsmc = int (nsmc)
      numc = int (numc)
      phelp = make_parallel_helper (parallel)
  
      if t.ndim != 1:
          raise ValueError ('`t` must be <= 1D')
  
      if x.shape != t.shape:
          raise ValueError ('`t` and `x` arguments must be the same size')
  
      if u.shape != t.shape:
          raise ValueError ('`t` and `u` arguments must be the same size')
  
      if periods.ndim != 1:
          raise ValueError ('`periods` must be <= 1D')
  
      if nbin < 2:
          raise ValueError ('`nbin` must be at least 2')
  
      if nshift < 1:
          raise ValueError ('`nshift` must be at least 1')
  
      if nsmc < 0:
          raise ValueError ('`nsmc` must be nonnegative')
  
      if numc < 0:
          raise ValueError ('`numc` must be nonnegative')
  
      # We can finally get started!
  
      if weights:
          wt = u
          u = wt ** -0.5
      else:
          wt = u ** -2
  
      v_all = weighted_variance (x, wt)
  
      with phelp.get_map () as map:
          get_thetas = lambda args: np.asarray (map (_map_one_theta, args))
          thetas = get_thetas ((t, x, wt, p, nbin, nshift, v_all)
                               for p in periods)
          imin = thetas.argmin ()
          pmin = periods[imin]
  
          # Now do the Monte Carlo jacknifing so that the caller can have some idea
          # as to the significance of the minimal value of `thetas`.
  
          mc_thetas = np.empty (periods.shape)
          mc_tmins = np.empty (nsmc)
  
          for i in range (nsmc):
              shuf = np.random.permutation (x.size)
              # Note that what we do here is very MapReduce-y. I'm not aware of
              # an easy way to implement this computation in that model, though.
              mc_thetas = get_thetas ((t, x[shuf], wt[shuf], p, nbin, nshift, v_all)
                                      for p in periods)
              mc_tmins[i] = mc_thetas.min ()
  
          mc_tmins.sort ()
          mc_pvalue = mc_tmins.searchsorted (thetas[imin]) / nsmc
  
          # Now add noise to assess the uncertainty of the period.
  
          mc_pmins = np.empty (numc)
  
          for i in range (numc):
              noised = np.random.normal (x, u)
              mc_thetas = get_thetas ((t, noised, wt, p, nbin, nshift, v_all)
                                      for p in periods)
              mc_pmins[i] = periods[mc_thetas.argmin ()]
  
          mc_pmins.sort ()
          mc_puncert = mc_pmins.std ()
  
      # All done.
  
      return PDMResult (thetas=thetas, imin=imin, pmin=pmin, mc_tmins=mc_tmins,
                        mc_pvalue=mc_pvalue, mc_pmins=mc_pmins,
                        mc_puncert=mc_puncert)