def surface_deriv(self, param):
        """Return the surface derivative at ``param``.

        The derivative at parameter ``phi`` is given by ::

            deriv = R * radius * (-sin(phi), -cos(phi))

        where R is a rotation matrix.

        Parameters
        ----------
        param : float or `array-like`
            Parameter value(s) at which to evaluate.

        Returns
        -------
        deriv : `numpy.ndarray`
            Array representing the derivative vector(s) at ``param``.
            If ``param`` is a single parameter, the returned array has
            shape ``(2,)``, otherwise ``param.shape + (2,)``.

        See Also
        --------
        surface

        Examples
        --------
        The method works with a single parameter, resulting in a single
        vector:

        >>> part = odl.uniform_partition(-np.pi / 2, np.pi / 2, 10)
        >>> det = CircularDetector(part, axis=[1, 0], radius=2)
        >>> det.surface_deriv(0)
        array([ 2.,  0.])

        It is also vectorized, i.e., it can be called with multiple
        parameters at once (or an n-dimensional array of parameters):

        >>> np.round(det.surface_deriv([-np.pi / 2, 0, np.pi / 2]), 10)
        array([[ 0.,  2.],
               [ 2.,  0.],
               [ 0., -2.]])

        >>> det.surface_deriv(np.zeros((4, 5))).shape
        (4, 5, 2)
        """
        squeeze_out = (np.shape(param) == ())
        param = np.array(param, dtype=float, copy=False, ndmin=1)
        if self.check_bounds and not is_inside_bounds(param, self.params):
            raise ValueError('`param` {} not in the valid range '
                             '{}'.format(param, self.params))

        deriv = np.empty(param.shape + (2,))
        deriv[..., 0] = -np.sin(param)
        deriv[..., 1] = -np.cos(param)
        deriv *= self.radius
        deriv = np.matmul(deriv, np.transpose(self.rotation_matrix))

        if squeeze_out:
            deriv = deriv.squeeze()

        return deriv