[python-dtcwt] 113/497: a first stab at keypoint detectors
Ghislain Vaillant
ghisvail-guest at moszumanska.debian.org
Tue Jul 21 18:05:55 UTC 2015
This is an automated email from the git hooks/post-receive script.
ghisvail-guest pushed a commit to branch debian/sid
in repository python-dtcwt.
commit 49a03a74e50722e76b550c38c193a3e510172128
Author: Rich Wareham <rjw57 at cam.ac.uk>
Date: Tue Aug 20 16:33:13 2013 +0100
a first stab at keypoint detectors
---
docs/reference.rst | 6 ++
dtcwt/keypoint.py | 245 +++++++++++++++++++++++++++++++++++++++++++++++++++++
2 files changed, 251 insertions(+)
diff --git a/docs/reference.rst b/docs/reference.rst
index f4f0d9b..9022605 100644
--- a/docs/reference.rst
+++ b/docs/reference.rst
@@ -7,6 +7,12 @@ Computing the DT-CWT
.. automodule:: dtcwt
:members:
+Keypoint analysis
+`````````````````
+
+.. automodule:: dtcwt.keypoint
+ :members:
+
Image sampling
``````````````
diff --git a/dtcwt/keypoint.py b/dtcwt/keypoint.py
new file mode 100644
index 0000000..91437fe
--- /dev/null
+++ b/dtcwt/keypoint.py
@@ -0,0 +1,245 @@
+import numpy as np
+
+from dtcwt.sampling import upsample_highpass, upsample
+
+__all__ = [ 'find_keypoints' ]
+
+def find_keypoints(highpass_subbands, method=None, alpha=1.0, beta=0.4,
+ threshold=None, max_points=None,
+ upsample_keypoint_energy=None, upsample_subbands=None,
+ refine_positions=True, skip_levels=1):
+ """
+ :param highpass_subbands: (NxMx6) matrix of highpass subband images
+ :param method: *(optional)* string specifying which keypoint energy method to use
+ :param alpha: *(optional)* scale parameter for ``'fauqueur'`` method
+ :param beta: *(optional)* shape parameter for ``'fauqueur'`` method
+ :param threshold: *(optional)* minimum keypoint energy of returned keypoints
+ :param max_points: *(optional)* maximum number of keypoints to return
+ :param upsample_keypoint_energy: is non-None, a string specifying a method used to upscale the keypoint energy map before finding keypoints
+ :param upsample_subands: is non-None, a string specifying a method used to upscale the subband image before finding keypoints
+ :param refine_positions: *(optional)* should the keypoint positions be refined to sub-pixel accuracy
+ :param skip_levels: *(optional)* number of levels of the transform to ignore before looking for keypoints
+
+ :returns: (Px4) array of P keypoints in image co-ordinates
+
+ The rows of the returned keypoint array give the x co-ordinate, y
+ co-ordinate, scale and keypoint energy. The rows are sorted in order of
+ decreasing keypoint energy.
+
+ If *refine_positions* is ``True`` then the positions (and energy) of the
+ keypoints will be refined to sub-pixel accuracy by fitting a quadratic
+ patch. If *refine_positions* is ``False`` then the keypoint locations will
+ be those corresponding directly to pixel-wise maxima of the subband images.
+
+ The *max_points* and *threshold* parameters are cumulative: if both are
+ specified then the *max_points* greatest energy keypoints with energy
+ greater than *threshold* will be returned.
+
+ Usually the keypoint energies returned from the finest scale level are
+ dominated by noise and so one usually wants to specify *skip_levels* to be
+ 1 or 2. If *skip_levels* is 0 then all levels will be used to compute
+ keypoint energy.
+
+ The *upsample_subbands* and *upsample_keypoint_energy* parameters are used
+ to control whether the individual subband coefficients and/org the keypoint
+ energy map are upscaled by 2 before finding keypoints. If these parameters
+ are None then no corresponding upscaling is performed. If non-None they
+ specify the upscale method as outlined in
+ :py:func:`dtcwt.sampling.upsample`.
+
+ If *method* is ``None``, the default ``'fauqueur'`` method is used.
+
+ =========== ======================================= ======================
+ Name Description Parameters used
+ =========== ======================================= ======================
+ fauqueur Geometric mean of absolute values[1] *alpha*, *beta*
+ bendale Minimum absolute value[2] none
+ kingsbury Cross-product of orthogonal subbands none
+ gale Gradient of subbands none
+ =========== ======================================= ======================
+
+ [1] Julien Fauqueur, Nick Kingsbury, and Ryan Anderson. *Multiscale
+ Keypoint Detection using the Dual-Tree Complex Wavelet Transform*. 2006
+ International Conference on Image Processing, pages 1625-1628, October
+ 2006. ISSN 1522-4880. doi: 10.1109/ICIP.2006.312656.
+ http://ieeexplore.ieee.org/lpdocs/epic03/wrapper.htm?arnumber=4106857.
+
+ [2] Pashmina Bendale, Bill Triggs, and Nick Kingsbury. *Multiscale Keypoint
+ Analysis based on Complex Wavelets*. In British Machine Vision Con-ference
+ (BMVC), 2010.
+ http://www-sigproc.eng.cam.ac.uk/~pb397/publications/BTK_BMVC_2010_abstract.pdf.
+ """
+ # Set default method
+ if method is None:
+ method = 'fauqueur'
+
+ # Skip levels
+ highpass_subbands = highpass_subbands[skip_levels:]
+
+ # Compute contribution to scale from upsampling
+ upsample_scale = 1
+ if upsample_subbands is not None:
+ upsample_scale <<= 1
+ if upsample_keypoint_energy is not None:
+ upsample_scale <<= 1
+
+ # Find keypoint energy map for each level
+ kp_energies = []
+ for subband in highpass_subbands:
+ if upsample_subbands is not None:
+ subband = upsample_highpass(subband, upsample_subbands)
+
+ if method == 'fauqueur':
+ kp_energies.append(_keypoint_energy_fauqueur(subband, alpha, beta))
+ elif method == 'bendale':
+ kp_energies.append(_keypoint_energy_bendale(subband))
+ elif method == 'kingsbury':
+ kp_energies.append(_keypoint_energy_kingsbury(subband))
+ elif method == 'gale':
+ kp_energies.append(_keypoint_energy_gale(subband))
+ else:
+ raise ValueError('Unknown method: {0}'.format(method))
+
+ if upsample_keypoint_energy is not None:
+ kp_energies[-1] = upsample(kp_energies[-1], upsample_keypoint_energy)
+
+ # Find keypoints for each level
+ kps = None
+ for level_idx, kp_energy in enumerate(kp_energies):
+ kp_scale = 2**(level_idx+1+skip_levels) / float(upsample_scale)
+ kp_rows, kp_cols, kp_energies = _kp_energy_maxima(kp_energy, threshold=threshold, refine=refine_positions)
+
+ # Scaling is a bit non-trivial. If the subband has pixel coords in range {0, .., M-1} then it has extent
+ # (-0.5, M-0.5]. If we need to scale the pixel size by kp_scale then the final image will have extent
+ # (-0.5, kp_scale*M-0.5]. So we need a linear function which maps -0.5 -> -0.5 and M-0.5 -> kp_scale*M-0.5
+ # such a function is x -> kp_scale * (x+0.5) - 0.5
+
+ level_kps = np.array((
+ (kp_cols+0.5)*kp_scale-0.5, (kp_rows+0.5)*kp_scale-0.5,
+ kp_scale*np.ones(kp_cols.shape[0]), kp_energies)).T
+
+ if kps is None:
+ kps = level_kps
+ else:
+ kps = np.vstack((kps, level_kps))
+
+ # Sort keypoints
+ sorted_indices = np.argsort(kps[:, 3])
+ kps = kps[sorted_indices[::-1],:]
+
+ # Truncate if necessary
+ if max_points is not None:
+ kps = kps[:max_points]
+
+ # Return keypoints
+ return kps
+
+def _keypoint_energy_fauqueur(subband, alpha, beta):
+ return alpha * np.power(np.maximum(0, np.product(np.abs(subband), axis=2)), beta)
+
+def _keypoint_energy_bendale(subband):
+ return np.min(np.abs(subband), axis=2)
+
+def _keypoint_energy_kingsbury(subband):
+ raise NotImplementedError('not implemented yet')
+
+def _keypoint_energy_gale(subband):
+ raise NotImplementedError('not implemented yet')
+
+def _nullspace(A, atol=1e-13, rtol=0):
+ """Compute an approximate basis for the nullspace of A.
+
+ The algorithm used by this function is based on the singular value
+ decomposition of `A`.
+
+ Parameters
+ ----------
+ A : ndarray
+ A should be at most 2-D. A 1-D array with length k will be treated
+ as a 2-D with shape (1, k)
+ atol : float
+ The absolute tolerance for a zero singular value. Singular values
+ smaller than `atol` are considered to be zero.
+ rtol : float
+ The relative tolerance. Singular values less than rtol*smax are
+ considered to be zero, where smax is the largest singular value.
+
+ If both `atol` and `rtol` are positive, the combined tolerance is the
+ maximum of the two; that is::
+ tol = max(atol, rtol * smax)
+ Singular values smaller than `tol` are considered to be zero.
+
+ Return value
+ ------------
+ ns : ndarray
+ If `A` is an array with shape (m, k), then `ns` will be an array
+ with shape (k, n), where n is the estimated dimension of the
+ nullspace of `A`. The columns of `ns` are a basis for the
+ nullspace; each element in numpy.dot(A, ns) will be approximately
+ zero.
+ """
+
+ A = np.atleast_2d(A)
+ u, s, vh = np.linalg.svd(A)
+ tol = max(atol, rtol * s[0])
+ nnz = (s >= tol).sum()
+ ns = vh[nnz:].conj().T
+ return ns
+
+def _kp_energy_maxima(X, threshold=None, refine=True):
+ # If no threshold is provided, choose one which all keypoints will pass
+ if threshold is None:
+ threshold = X.min() - 1
+
+ # Compute local maximum image
+ maxima = np.ones_like(X) * threshold
+ for dx in (-1,0,1):
+ for dy in (-1,0,1):
+ maxima[1:-2,1:-2] = np.maximum(maxima[1:-2,1:-2], X[1+dy:X.shape[0]-2+dy, 1+dx:X.shape[1]-2+dx])
+
+ # This will be used to store the values of local maxima
+ lm_values = []
+
+ # This will be used to store the refined positions
+ lm_refined_rows, lm_refined_cols = [], []
+
+ # Find local maxima
+ lm_rows, lm_cols = np.nonzero(maxima == X)
+
+ if refine:
+ # Taylor series of I(x) around x_0 is I(x) ~= I(x_o) + dI/dx x + dI^2/dx^2 x^2 + ...
+ # maximum is at differential is zero or:
+ # 0 = dI/dx + 2 dI^2/dx^2 x => x = -dI/dx * (2*dI^2/dx^2)^-1
+
+ # Form the various gradient images for X
+ dXdy, dXdx = np.gradient(X)
+ dX2dxdy, dX2dx2 = np.gradient(dXdx)
+ dX2dy2, dX2dydx = np.gradient(dXdy)
+
+ a_im = np.dstack((
+ dX2dx2, dX2dy2, dX2dxdy, dXdx, dXdy, X,
+ ))
+
+ # Calculate a vectors for each neighbourhood
+ for r, c in zip(lm_rows, lm_cols):
+ if refine:
+ a = a_im[r,c,:]
+ A = np.array(((2*a[0], a[2], a[3]), (a[2], 2*a[1], a[4])))
+ v = _nullspace(A)[:,0]
+ v /= v[2]
+
+ # Only accept fittings where maximum is within half of a pixel of
+ # the maximum pixel's centre.
+ if np.any(np.abs(v[:2]) > 0.5):
+ continue
+
+ x, y = v[:2]
+ lm_values.append(a[0]*x*x + a[1]*y*y + a[2]*x*y + a[3]*x + a[4]*y + a[5])
+ else:
+ x, y = 0, 0
+ lm_values.append(X[r,c])
+
+ lm_refined_rows.append(r+y)
+ lm_refined_cols.append(c+x)
+
+ return np.array(lm_refined_rows), np.array(lm_refined_cols), np.array(lm_values)
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-science/packages/python-dtcwt.git
More information about the debian-science-commits
mailing list