mirror of
https://github.com/wassname/scikit-image.git
synced 2026-08-01 12:50:48 +08:00
Merge radon/fast_homography.
This commit is contained in:
Stefan van der Walt
+5
-1
@@ -44,7 +44,8 @@
|
||||
Incorporating CellProfiler's Sobel edge detector, build and bug fixes.
|
||||
|
||||
- Emmanuelle Guillart
|
||||
Total variation noise filtering
|
||||
Total variation noise filtering, integration of CellProfiler's
|
||||
mathematical morphology tools.
|
||||
|
||||
- Maël Primet
|
||||
Total variation noise filtering
|
||||
@@ -60,3 +61,6 @@
|
||||
|
||||
- Kyle Mandli
|
||||
CSV to ReST code for feature comparison table.
|
||||
|
||||
- The Scikit Learn team
|
||||
From whom we borrowed the example generation tools.
|
||||
|
||||
@@ -0,0 +1,6 @@
|
||||
|
||||
General examples
|
||||
-------------------
|
||||
|
||||
General-purpose and introductory examples for the scikit.
|
||||
|
||||
@@ -0,0 +1,52 @@
|
||||
"""
|
||||
====================================================
|
||||
Denoising the picture of Lena using total variation
|
||||
====================================================
|
||||
|
||||
In this example, we denoise a noisy version of the picture of Lena
|
||||
using the total variation denoising filter. The result of this filter
|
||||
is an image that has a minimal total variation norm, while being as
|
||||
close to the initial image as possible. The total variation is the L1
|
||||
norm of the gradient of the image, and minimizing the total variation
|
||||
typically produces "posterized" images with flat domains separated by
|
||||
sharp edges.
|
||||
|
||||
It is possible to change the degree of posterization by controlling
|
||||
the tradeoff between denoising and faithfulness to the original image.
|
||||
|
||||
"""
|
||||
|
||||
import numpy as np
|
||||
import matplotlib.pyplot as plt
|
||||
|
||||
from scikits.image import data
|
||||
from scikits.image.filter import tv_denoise
|
||||
|
||||
l = data.lena()
|
||||
l = l[230:290, 220:320]
|
||||
|
||||
noisy = l + 0.4*l.std()*np.random.random(l.shape)
|
||||
|
||||
tv_denoised = tv_denoise(noisy, weight=10)
|
||||
|
||||
|
||||
plt.figure(figsize=(12,2.8))
|
||||
|
||||
plt.subplot(131)
|
||||
plt.imshow(noisy, cmap=plt.cm.gray, vmin=40, vmax=220)
|
||||
plt.axis('off')
|
||||
plt.title('noisy', fontsize=20)
|
||||
plt.subplot(132)
|
||||
plt.imshow(tv_denoised, cmap=plt.cm.gray, vmin=40, vmax=220)
|
||||
plt.axis('off')
|
||||
plt.title('TV denoising', fontsize=20)
|
||||
|
||||
tv_denoised = tv_denoise(noisy, weight=50)
|
||||
plt.subplot(133)
|
||||
plt.imshow(tv_denoised, cmap=plt.cm.gray, vmin=40, vmax=220)
|
||||
plt.axis('off')
|
||||
plt.title('(more) TV denoising', fontsize=20)
|
||||
|
||||
plt.subplots_adjust(wspace=0.02, hspace=0.02, top=0.9, bottom=0, left=0,
|
||||
right=1)
|
||||
plt.show()
|
||||
@@ -66,3 +66,46 @@ products or services of Licensee, or any third party.
|
||||
8. By copying, installing or otherwise using matplotlib 0.98.3,
|
||||
Licensee agrees to be bound by the terms and conditions of this
|
||||
License Agreement.
|
||||
|
||||
The file
|
||||
|
||||
- gen_rst.py
|
||||
|
||||
was taken from the scikit-learn (http://scikit-learn.sourceforge.net), which
|
||||
has the following license:
|
||||
|
||||
New BSD License
|
||||
|
||||
Copyright (c) 2007 - 2011 The scikit-learn developers.
|
||||
All rights reserved.
|
||||
|
||||
|
||||
Redistribution and use in source and binary forms, with or without
|
||||
modification, are permitted provided that the following conditions are met:
|
||||
|
||||
a. Redistributions of source code must retain the above copyright notice,
|
||||
this list of conditions and the following disclaimer.
|
||||
b. Redistributions in binary form must reproduce the above copyright
|
||||
notice, this list of conditions and the following disclaimer in the
|
||||
documentation and/or other materials provided with the distribution.
|
||||
c. Neither the name of the Scikit-learn Developers nor the names of
|
||||
its contributors may be used to endorse or promote products
|
||||
derived from this software without specific prior written
|
||||
permission.
|
||||
|
||||
|
||||
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
|
||||
AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
|
||||
IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
|
||||
ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE FOR
|
||||
ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
|
||||
DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
|
||||
SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
|
||||
CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
|
||||
LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
|
||||
OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
|
||||
DAMAGE.
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
@@ -0,0 +1,305 @@
|
||||
"""
|
||||
Example generation for the scikit image.
|
||||
|
||||
Generate the rst files for the examples by iterating over the python
|
||||
example files.
|
||||
|
||||
Files that generate images should start with 'plot'.
|
||||
|
||||
This code was taken from the scikit-learn.
|
||||
|
||||
"""
|
||||
import os
|
||||
import shutil
|
||||
import traceback
|
||||
import glob
|
||||
|
||||
import matplotlib
|
||||
matplotlib.use('Agg')
|
||||
|
||||
import token, tokenize
|
||||
|
||||
rst_template = """
|
||||
|
||||
.. _example_%(short_fname)s:
|
||||
|
||||
%(docstring)s
|
||||
|
||||
**Python source code:** :download:`%(fname)s <%(fname)s>`
|
||||
|
||||
.. literalinclude:: %(fname)s
|
||||
:lines: %(end_row)s-
|
||||
"""
|
||||
|
||||
plot_rst_template = """
|
||||
|
||||
.. _example_%(short_fname)s:
|
||||
|
||||
%(docstring)s
|
||||
|
||||
%(image_list)s
|
||||
|
||||
**Python source code:** :download:`%(fname)s <%(fname)s>`
|
||||
|
||||
.. literalinclude:: %(fname)s
|
||||
:lines: %(end_row)s-
|
||||
"""
|
||||
|
||||
# The following strings are used when we have several pictures: we use
|
||||
# an html div tag that our CSS uses to turn the lists into horizontal
|
||||
# lists.
|
||||
HLIST_HEADER = """
|
||||
.. rst-class:: horizontal
|
||||
|
||||
"""
|
||||
|
||||
HLIST_IMAGE_TEMPLATE = """
|
||||
*
|
||||
|
||||
.. image:: images/%s
|
||||
:scale: 50
|
||||
"""
|
||||
|
||||
SINGLE_IMAGE = """
|
||||
.. image:: images/%s
|
||||
:align: center
|
||||
"""
|
||||
|
||||
def extract_docstring(filename):
|
||||
""" Extract a module-level docstring, if any
|
||||
"""
|
||||
lines = file(filename).readlines()
|
||||
start_row = 0
|
||||
if lines[0].startswith('#!'):
|
||||
lines.pop(0)
|
||||
start_row = 1
|
||||
|
||||
docstring = ''
|
||||
first_par = ''
|
||||
tokens = tokenize.generate_tokens(lines.__iter__().next)
|
||||
for tok_type, tok_content, _, (erow, _), _ in tokens:
|
||||
tok_type = token.tok_name[tok_type]
|
||||
if tok_type in ('NEWLINE', 'COMMENT', 'NL', 'INDENT', 'DEDENT'):
|
||||
continue
|
||||
elif tok_type == 'STRING':
|
||||
docstring = eval(tok_content)
|
||||
# If the docstring is formatted with several paragraphs, extract
|
||||
# the first one:
|
||||
paragraphs = '\n'.join(line.rstrip()
|
||||
for line in docstring.split('\n')).split('\n\n')
|
||||
if len(paragraphs) > 0:
|
||||
first_par = paragraphs[0]
|
||||
break
|
||||
return docstring, first_par, erow+1+start_row
|
||||
|
||||
|
||||
def generate_example_rst(app):
|
||||
""" Generate the list of examples, as well as the contents of
|
||||
examples.
|
||||
"""
|
||||
root_dir = os.path.join(app.builder.srcdir, 'auto_examples')
|
||||
example_dir = os.path.abspath(app.builder.srcdir + '/../' + 'examples')
|
||||
try:
|
||||
plot_gallery = eval(app.builder.config.plot_gallery)
|
||||
except TypeError:
|
||||
plot_gallery = bool(app.builder.config.plot_gallery)
|
||||
if not os.path.exists(example_dir):
|
||||
os.makedirs(example_dir)
|
||||
if not os.path.exists(root_dir):
|
||||
os.makedirs(root_dir)
|
||||
|
||||
# we create an index.rst with all examples
|
||||
fhindex = file(os.path.join(root_dir, 'index.txt'), 'w')
|
||||
fhindex.write("""\
|
||||
|
||||
.. raw:: html
|
||||
|
||||
<style type="text/css">
|
||||
.figure {
|
||||
float: left;
|
||||
margin: 10px;
|
||||
width: auto;
|
||||
height: 200px;
|
||||
width: 180px;
|
||||
}
|
||||
|
||||
.figure img {
|
||||
display: inline;
|
||||
}
|
||||
|
||||
.figure .caption {
|
||||
width: 170px;
|
||||
text-align: center !important;
|
||||
}
|
||||
</style>
|
||||
|
||||
Examples
|
||||
========
|
||||
|
||||
.. _examples-index:
|
||||
""")
|
||||
# Here we don't use an os.walk, but we recurse only twice: flat is
|
||||
# better than nested.
|
||||
generate_dir_rst('.', fhindex, example_dir, root_dir, plot_gallery)
|
||||
for dir in sorted(os.listdir(example_dir)):
|
||||
if os.path.isdir(os.path.join(example_dir, dir)):
|
||||
generate_dir_rst(dir, fhindex, example_dir, root_dir, plot_gallery)
|
||||
fhindex.flush()
|
||||
|
||||
|
||||
def generate_dir_rst(dir, fhindex, example_dir, root_dir, plot_gallery):
|
||||
""" Generate the rst file for an example directory.
|
||||
"""
|
||||
if not dir == '.':
|
||||
target_dir = os.path.join(root_dir, dir)
|
||||
src_dir = os.path.join(example_dir, dir)
|
||||
else:
|
||||
target_dir = root_dir
|
||||
src_dir = example_dir
|
||||
if not os.path.exists(os.path.join(src_dir, 'README.txt')):
|
||||
print 80*'_'
|
||||
print ('Example directory %s does not have a README.txt file'
|
||||
% src_dir)
|
||||
print 'Skipping this directory'
|
||||
print 80*'_'
|
||||
return
|
||||
fhindex.write("""
|
||||
|
||||
|
||||
%s
|
||||
|
||||
|
||||
""" % file(os.path.join(src_dir, 'README.txt')).read())
|
||||
if not os.path.exists(target_dir):
|
||||
os.makedirs(target_dir)
|
||||
|
||||
def sort_key(a):
|
||||
# put last elements without a plot
|
||||
if not a.startswith('plot') and a.endswith('.py'):
|
||||
return 'zz' + a
|
||||
return a
|
||||
for fname in sorted(os.listdir(src_dir), key=sort_key):
|
||||
if fname.endswith('py'):
|
||||
generate_file_rst(fname, target_dir, src_dir, plot_gallery)
|
||||
thumb = os.path.join(dir, 'images', 'thumb', fname[:-3] + '.png')
|
||||
link_name = os.path.join(dir, fname).replace(os.path.sep, '_')
|
||||
fhindex.write('.. figure:: %s\n' % thumb)
|
||||
if link_name.startswith('._'):
|
||||
link_name = link_name[2:]
|
||||
if dir != '.':
|
||||
fhindex.write(' :target: ./%s/%s.html\n\n' % (dir, fname[:-3]))
|
||||
else:
|
||||
fhindex.write(' :target: ./%s.html\n\n' % link_name[:-3])
|
||||
fhindex.write(' :ref:`example_%s`\n\n' % link_name)
|
||||
fhindex.write("""
|
||||
.. raw:: html
|
||||
|
||||
<div style="clear: both"></div>
|
||||
""") # clear at the end of the section
|
||||
|
||||
|
||||
def generate_file_rst(fname, target_dir, src_dir, plot_gallery):
|
||||
""" Generate the rst file for a given example.
|
||||
"""
|
||||
base_image_name = os.path.splitext(fname)[0]
|
||||
image_fname = '%s_%%s.png' % base_image_name
|
||||
|
||||
this_template = rst_template
|
||||
last_dir = os.path.split(src_dir)[-1]
|
||||
# to avoid leading . in file names, and wrong names in links
|
||||
if last_dir == '.' or last_dir == 'examples':
|
||||
last_dir = ''
|
||||
else:
|
||||
last_dir += '_'
|
||||
short_fname = last_dir + fname
|
||||
src_file = os.path.join(src_dir, fname)
|
||||
example_file = os.path.join(target_dir, fname)
|
||||
shutil.copyfile(src_file, example_file)
|
||||
|
||||
# The following is a list containing all the figure names
|
||||
figure_list = []
|
||||
|
||||
image_dir = os.path.join(target_dir, 'images')
|
||||
thumb_dir = os.path.join(image_dir, 'thumb')
|
||||
if not os.path.exists(image_dir):
|
||||
os.makedirs(image_dir)
|
||||
if not os.path.exists(thumb_dir):
|
||||
os.makedirs(thumb_dir)
|
||||
image_path = os.path.join(image_dir, image_fname)
|
||||
thumb_file = os.path.join(thumb_dir, fname[:-3] + '.png')
|
||||
if plot_gallery and fname.startswith('plot'):
|
||||
# generate the plot as png image if file name
|
||||
# starts with plot and if it is more recent than an
|
||||
# existing image.
|
||||
first_image_file = image_path % 1
|
||||
|
||||
if (not os.path.exists(first_image_file) or
|
||||
os.stat(first_image_file).st_mtime <=
|
||||
os.stat(src_file).st_mtime):
|
||||
# We need to execute the code
|
||||
print 'plotting %s' % fname
|
||||
import matplotlib.pyplot as plt
|
||||
plt.close('all')
|
||||
cwd = os.getcwd()
|
||||
try:
|
||||
# First CD in the original example dir, so that any file created
|
||||
# by the example get created in this directory
|
||||
os.chdir(os.path.dirname(src_file))
|
||||
execfile(os.path.basename(src_file), {'pl' : plt})
|
||||
os.chdir(cwd)
|
||||
|
||||
# In order to save every figure we have two solutions :
|
||||
# * iterate from 1 to infinity and call plt.fignum_exists(n)
|
||||
# (this requires the figures to be numbered
|
||||
# incrementally: 1, 2, 3 and not 1, 2, 5)
|
||||
# * iterate over [fig_mngr.num for fig_mngr in
|
||||
# matplotlib._pylab_helpers.Gcf.get_all_fig_managers()]
|
||||
for fig_num in (fig_mngr.num for fig_mngr in
|
||||
matplotlib._pylab_helpers.Gcf.get_all_fig_managers()):
|
||||
# Set the fig_num figure as the current figure as we can't
|
||||
# save a figure that's not the current figure.
|
||||
plt.figure(fig_num)
|
||||
plt.savefig(image_path % fig_num)
|
||||
figure_list.append(image_fname % fig_num)
|
||||
except:
|
||||
print 80*'_'
|
||||
print '%s is not compiling:' % fname
|
||||
traceback.print_exc()
|
||||
print 80*'_'
|
||||
finally:
|
||||
os.chdir(cwd)
|
||||
else:
|
||||
figure_list = [f[len(image_dir):]
|
||||
for f in glob.glob(image_path % '[1-9]')]
|
||||
#for f in glob.glob(image_path % '*')]
|
||||
|
||||
# generate thumb file
|
||||
this_template = plot_rst_template
|
||||
from matplotlib import image
|
||||
if os.path.exists(first_image_file):
|
||||
image.thumbnail(first_image_file, thumb_file, 0.2)
|
||||
|
||||
if not os.path.exists(thumb_file):
|
||||
# create something not to replace the thumbnail
|
||||
shutil.copy('images/blank_image.png', thumb_file)
|
||||
|
||||
docstring, short_desc, end_row = extract_docstring(example_file)
|
||||
|
||||
# Depending on whether we have one or more figures, we're using a
|
||||
# horizontal list or a single rst call to 'image'.
|
||||
if len(figure_list) == 1:
|
||||
figure_name = figure_list[0]
|
||||
image_list = SINGLE_IMAGE % figure_name.lstrip('/')
|
||||
else:
|
||||
image_list = HLIST_HEADER
|
||||
for figure_name in figure_list:
|
||||
image_list += HLIST_IMAGE_TEMPLATE % figure_name.lstrip('/')
|
||||
|
||||
f = open(os.path.join(target_dir, fname[:-2] + 'txt'),'w')
|
||||
f.write(this_template % locals())
|
||||
f.flush()
|
||||
|
||||
|
||||
def setup(app):
|
||||
app.connect('builder-inited', generate_example_rst)
|
||||
app.add_config_value('plot_gallery', True, 'html')
|
||||
Binary file not shown.
|
After Width: | Height: | Size: 1.5 KiB |
+7
-1
@@ -23,11 +23,17 @@ sys.path.append(os.path.join(curpath, '..', 'ext'))
|
||||
|
||||
# -- General configuration -----------------------------------------------------
|
||||
|
||||
try:
|
||||
import gen_rst
|
||||
except:
|
||||
pass
|
||||
|
||||
|
||||
# Add any Sphinx extension module names here, as strings. They can be extensions
|
||||
# coming with Sphinx (named 'sphinx.ext.*') or your custom ones.
|
||||
extensions = ['sphinx.ext.autodoc', 'sphinx.ext.pngmath', 'numpydoc',
|
||||
'sphinx.ext.autosummary', 'sphinx.ext.inheritance_diagram',
|
||||
'plot_directive']
|
||||
'plot_directive', 'gen_rst']
|
||||
|
||||
# Add any paths that contain templates here, relative to this directory.
|
||||
templates_path = ['_templates']
|
||||
|
||||
Regular → Executable
+5
-5
@@ -40,16 +40,16 @@ def calculate_coverage(reader):
|
||||
total_items += 1
|
||||
if ":done:" in row[2] or ":done:" in row[3]:
|
||||
done_items += 1
|
||||
elif ":partial:" in row[2] or ":partial:" in row[3]:
|
||||
if ":partial:" in row[2] or ":partial:" in row[3]:
|
||||
partial_items += 1
|
||||
elif ":na:" in row[2] or ":na:" in row[3]:
|
||||
if ":na:" in row[2] or ":na:" in row[3]:
|
||||
na_items += 1
|
||||
|
||||
counts = (done_items,
|
||||
partial_items,
|
||||
total_items - (partial_items + done_items + na_items),
|
||||
na_items)
|
||||
|
||||
|
||||
return list(i / total_items for i in counts)
|
||||
|
||||
def read_table_titles(reader):
|
||||
@@ -235,7 +235,7 @@ Coverage Bar
|
||||
|
||||
for item, style in enumerate(('done-bar', 'partial-bar',
|
||||
'missing-bar', 'na-bar')):
|
||||
stream.write('<td width="%s" class="%s">XX</td>' % \
|
||||
stream.write('<td width="%s%%" class="%s"> </td>' % \
|
||||
(item_counts[item] * 100, style))
|
||||
|
||||
stream.write("</tr></table>\n\n")
|
||||
@@ -269,6 +269,6 @@ if __name__ == "__main__":
|
||||
generate_page(csv_path,output)
|
||||
output.close()
|
||||
|
||||
print "Generated %s from %s." % (output_path,csv_path)
|
||||
print("Generated %s from %s." % (output_path,csv_path))
|
||||
|
||||
|
||||
|
||||
@@ -39,6 +39,11 @@ Sections
|
||||
Conditions on the use and redistribution of this package.
|
||||
|
||||
|
||||
|
||||
`Examples <auto_examples/index.html>`_
|
||||
|
||||
Introductory examples
|
||||
|
||||
Indices and Tables
|
||||
==================
|
||||
|
||||
|
||||
@@ -0,0 +1,4 @@
|
||||
import scipy.misc as _m
|
||||
|
||||
lena = _m.lena
|
||||
|
||||
@@ -11,6 +11,7 @@ def configuration(parent_package='', top_path=None):
|
||||
config.add_subpackage('morphology')
|
||||
config.add_subpackage('filter')
|
||||
config.add_subpackage('transform')
|
||||
config.add_subpackage('data')
|
||||
|
||||
def add_test_directories(arg, dirname, fnames):
|
||||
if dirname.split(os.path.sep)[-1] == 'tests':
|
||||
|
||||
@@ -1,5 +1,6 @@
|
||||
from hough_transform import *
|
||||
from finite_radon_transform import *
|
||||
from radon_transform import *
|
||||
from project import *
|
||||
|
||||
from .hough_transform import *
|
||||
from .radon import *
|
||||
from .finite_radon_transform import *
|
||||
from .project import *
|
||||
from ._project import homography as fast_homography
|
||||
from .integral import *
|
||||
|
||||
@@ -0,0 +1,209 @@
|
||||
#cython: cdivison=True boundscheck=False
|
||||
|
||||
__all__ = ['homography']
|
||||
|
||||
cimport cython
|
||||
cimport numpy as np
|
||||
|
||||
import numpy as np
|
||||
import cython
|
||||
|
||||
from cython.operator import dereference
|
||||
|
||||
np.import_array()
|
||||
|
||||
cdef extern from "math.h":
|
||||
double floor(double)
|
||||
double fmod(double, double)
|
||||
|
||||
cdef double get_pixel(double *image, int rows, int cols,
|
||||
int r, int c, char mode, double cval=0):
|
||||
"""Get a pixel from the image, taking wrapping mode into consideration.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
image : *double
|
||||
Input image.
|
||||
rows, cols : int
|
||||
Dimensions of image.
|
||||
r, c : int
|
||||
Position at which to get the pixel.
|
||||
mode : {'C', 'W', 'M'}
|
||||
Wrapping mode. Constant, Wrap or Mirror.
|
||||
cval : double
|
||||
Constant value to use for mode constant.
|
||||
|
||||
"""
|
||||
if mode == 'C':
|
||||
if (r < 0) or (r > rows - 1) or (c < 0) or (c > cols - 1):
|
||||
return cval
|
||||
else:
|
||||
return image[r * cols + c]
|
||||
else:
|
||||
return image[coord_map(rows, r, mode) * cols +
|
||||
coord_map(cols, c, mode)]
|
||||
|
||||
cdef int coord_map(int dim, int coord, char mode):
|
||||
"""
|
||||
Wrap a coordinate, according to a given dimension and mode.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
dim : int
|
||||
Maximum coordinate.
|
||||
coord : int
|
||||
Coord provided by user. May be < 0 or > dim.
|
||||
mode : {'W', 'M'}
|
||||
Whether to wrap or mirror the coordinate if it
|
||||
falls outside [0, dim).
|
||||
|
||||
"""
|
||||
dim = dim - 1
|
||||
if mode == 'M': # mirror
|
||||
if (coord < 0):
|
||||
# How many times times does the coordinate wrap?
|
||||
if (<int>(-coord / dim) % 2 != 0):
|
||||
return dim - <int>(-coord % dim)
|
||||
else:
|
||||
return <int>(-coord % dim)
|
||||
elif (coord > dim):
|
||||
if (<int>(coord / dim) % 2 != 0):
|
||||
return <int>(dim - (coord % dim))
|
||||
else:
|
||||
return <int>(coord % dim)
|
||||
elif mode == 'W': # wrap
|
||||
if (coord < 0):
|
||||
return <int>(dim - (-coord % dim))
|
||||
elif (coord > dim):
|
||||
return <int>(coord % dim)
|
||||
|
||||
return coord
|
||||
|
||||
cdef tf(double x, double y, double* H, double *x_, double *y_):
|
||||
"""Apply a homography to a coordinate.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
x, y : double
|
||||
Input coordinate.
|
||||
H : (3,3) *double
|
||||
Transformation matrix.
|
||||
x_, y_ : *double
|
||||
Output coordinate.
|
||||
|
||||
"""
|
||||
cdef double xx, yy, zz
|
||||
|
||||
xx = H[0] * x + H[1] * y + H[2]
|
||||
yy = H[3] * x + H[4] * y + H[5]
|
||||
zz = H[6] * x + H[7] * y + H[8]
|
||||
|
||||
xx = xx / zz
|
||||
yy = yy / zz
|
||||
|
||||
x_[0] = xx
|
||||
y_[0] = yy
|
||||
|
||||
@cython.boundscheck(False)
|
||||
def homography(np.ndarray image, np.ndarray H, output_shape=None,
|
||||
mode='constant', double cval=0):
|
||||
"""
|
||||
Projective transformation (homography).
|
||||
|
||||
Perform a projective transformation (homography) of a
|
||||
floating point image, using bi-linear interpolation.
|
||||
|
||||
For each pixel, given its homogeneous coordinate :math:`\mathbf{x}
|
||||
= [x, y, 1]^T`, its target position is calculated by multiplying
|
||||
with the given matrix, :math:`H`, to give :math:`H \mathbf{x}`.
|
||||
E.g., to rotate by theta degrees clockwise, the matrix should be
|
||||
|
||||
::
|
||||
|
||||
[[cos(theta) -sin(theta) 0]
|
||||
[sin(theta) cos(theta) 0]
|
||||
[0 0 1]]
|
||||
|
||||
or, to translate x by 10 and y by 20,
|
||||
|
||||
::
|
||||
|
||||
[[1 0 10]
|
||||
[0 1 20]
|
||||
[0 0 1 ]].
|
||||
|
||||
Parameters
|
||||
----------
|
||||
image : 2-D array
|
||||
Input image.
|
||||
H : array of shape ``(3, 3)``
|
||||
Transformation matrix H that defines the homography.
|
||||
output_shape : tuple (rows, cols)
|
||||
Shape of the output image generated.
|
||||
order : int
|
||||
Order of splines used in interpolation.
|
||||
mode : {'constant', 'mirror', 'wrap'}
|
||||
How to handle values outside the image borders.
|
||||
cval : string
|
||||
Used in conjunction with mode 'C' (constant), the value
|
||||
outside the image boundaries.
|
||||
|
||||
"""
|
||||
|
||||
cdef np.ndarray[dtype=np.double_t, ndim=2] img = \
|
||||
np.asarray(image, dtype=np.double)
|
||||
cdef np.ndarray[dtype=np.double_t, ndim=2] M = \
|
||||
np.ascontiguousarray(np.linalg.inv(H))
|
||||
|
||||
if mode not in ('constant', 'wrap', 'mirror'):
|
||||
raise ValueError("Invalid mode specified. Please use "
|
||||
"`constant`, `wrap` or `mirror`.")
|
||||
if mode == 'constant':
|
||||
mode_c = ord('C')
|
||||
elif mode == 'wrap':
|
||||
mode_c = ord('W')
|
||||
elif mode == 'mirror':
|
||||
mode_c = ord('M')
|
||||
|
||||
cdef int out_r, out_c, columns, rows
|
||||
if output_shape is None:
|
||||
out_r = img.shape[0]
|
||||
out_c = img.shape[1]
|
||||
else:
|
||||
out_r = output_shape[0]
|
||||
out_c = output_shape[1]
|
||||
|
||||
rows = img.shape[0]
|
||||
columns = img.shape[1]
|
||||
|
||||
cdef np.ndarray[dtype=np.double_t, ndim=2] out = \
|
||||
np.zeros((out_r, out_c), dtype=np.double)
|
||||
|
||||
cdef int tfr, tfc, r_int, c_int
|
||||
cdef double y0, y1, y2, y3
|
||||
cdef double r, c, z, t, u
|
||||
|
||||
for tfr in range(out_r):
|
||||
for tfc in range(out_c):
|
||||
tf(tfc, tfr, <double*>M.data, &c, &r)
|
||||
r_int = <int>floor(r)
|
||||
c_int = <int>floor(c)
|
||||
|
||||
t = r - r_int
|
||||
u = c - c_int
|
||||
|
||||
y0 = get_pixel(<double*>img.data, rows, columns,
|
||||
r_int, c_int, mode_c)
|
||||
y1 = get_pixel(<double*>img.data, rows, columns,
|
||||
r_int + 1, c_int, mode_c)
|
||||
y2 = get_pixel(<double*>img.data, rows, columns,
|
||||
r_int + 1, c_int + 1, mode_c)
|
||||
y3 = get_pixel(<double*>img.data, rows, columns,
|
||||
r_int, c_int + 1, mode_c)
|
||||
|
||||
out[tfr, tfc] = \
|
||||
(1 - t) * (1 - u) * y0 + \
|
||||
t * (1 - u) * y1 + \
|
||||
t * u * y2 + (1 - t) * u * y3;
|
||||
|
||||
return out
|
||||
@@ -0,0 +1,60 @@
|
||||
def integral_image(x):
|
||||
"""Integral image / summed area table.
|
||||
|
||||
The integral image contains the sum of all elements above and to the
|
||||
left of it, i.e.:
|
||||
|
||||
.. math::
|
||||
|
||||
S[m, n] = \sum_{i \leq m} \sum_{j \leq n} X[i, j]
|
||||
|
||||
Parameters
|
||||
----------
|
||||
x : ndarray
|
||||
Input image.
|
||||
|
||||
Returns
|
||||
-------
|
||||
S : ndarray
|
||||
Integral image / summed area table.
|
||||
|
||||
References
|
||||
----------
|
||||
.. [1] F.C. Crow, "Summed-area tables for texture mapping,"
|
||||
ACM SIGGRAPH Computer Graphics, vol. 18, 1984, pp. 207-212.
|
||||
|
||||
"""
|
||||
return x.cumsum(1).cumsum(0)
|
||||
|
||||
def integrate(ii, r0, c0, r1, c1):
|
||||
"""Use an integral image to integrate over a given window.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
ii : ndarray
|
||||
Integral image.
|
||||
r0, c0 : int
|
||||
Top-left corner of block to be summed.
|
||||
r1, c1 : int
|
||||
Bottom-right corner of block to be summed.
|
||||
|
||||
Returns
|
||||
-------
|
||||
S : int
|
||||
Integral (sum) over the given window.
|
||||
|
||||
"""
|
||||
S = 0
|
||||
|
||||
S += ii[r1, c1]
|
||||
|
||||
if (r0 - 1 >= 0) and (c0 - 1 >= 0):
|
||||
S += ii[r0 - 1, c0 - 1]
|
||||
|
||||
if (r0 - 1 >= 0):
|
||||
S -= ii[r0 - 1, c1]
|
||||
|
||||
if (c0 - 1 >= 0):
|
||||
S -= ii[r1, c0 - 1]
|
||||
|
||||
return S
|
||||
@@ -15,10 +15,14 @@ def configuration(parent_package='', top_path=None):
|
||||
config.add_data_dir('tests')
|
||||
|
||||
cython(['_hough_transform.pyx'], working_path=base_path)
|
||||
cython(['_project.pyx'], working_path=base_path)
|
||||
|
||||
config.add_extension('_hough_transform', sources=['_hough_transform.c'],
|
||||
include_dirs=[get_numpy_include_dirs()])
|
||||
|
||||
config.add_extension('_project', sources=['_project.c'],
|
||||
include_dirs=[get_numpy_include_dirs()])
|
||||
|
||||
return config
|
||||
|
||||
if __name__ == '__main__':
|
||||
|
||||
@@ -0,0 +1,28 @@
|
||||
import numpy as np
|
||||
from numpy.testing import *
|
||||
|
||||
from scikits.image.transform import integral_image, integrate
|
||||
|
||||
x = (np.random.random((50, 50)) * 255).astype(np.uint8)
|
||||
s = integral_image(x)
|
||||
|
||||
def test_validity():
|
||||
y = np.arange(12).reshape((4,3))
|
||||
|
||||
y = (np.random.random((50, 50)) * 255).astype(np.uint8)
|
||||
assert_equal(integral_image(y)[-1, -1],
|
||||
y.sum())
|
||||
|
||||
def test_basic():
|
||||
assert_equal(x[12:24, 10:20].sum(), integrate(s, 12, 10, 23, 19))
|
||||
assert_equal(x[:20, :20].sum(), integrate(s, 0, 0, 19, 19))
|
||||
assert_equal(x[:20, 10:20].sum(), integrate(s, 0, 10, 19, 19))
|
||||
assert_equal(x[10:20, :20].sum(), integrate(s, 10, 0, 19, 19))
|
||||
|
||||
def test_single():
|
||||
assert_equal(x[0, 0], integrate(s, 0, 0, 0, 0))
|
||||
assert_equal(x[10, 10], integrate(s, 10, 10, 10, 10))
|
||||
|
||||
|
||||
if __name__ == '__main__':
|
||||
run_module_suite()
|
||||
@@ -1,7 +1,9 @@
|
||||
import numpy as np
|
||||
from numpy.testing import assert_array_almost_equal
|
||||
|
||||
from scikits.image.transform.project import _stackcopy, homography
|
||||
from scikits.image.transform.project import _stackcopy
|
||||
from scikits.image.transform import homography, fast_homography
|
||||
from scikits.image import data
|
||||
|
||||
def test_stackcopy():
|
||||
layers = 4
|
||||
@@ -19,3 +21,41 @@ def test_homography():
|
||||
[0, 0, 1]])
|
||||
x90 = homography(x, M, order=1)
|
||||
assert_array_almost_equal(x90, np.rot90(x))
|
||||
|
||||
def test_fast_homography():
|
||||
img = data.lena()
|
||||
img = img[:, :100]
|
||||
|
||||
theta = np.deg2rad(30)
|
||||
scale = 0.5
|
||||
tx, ty = 50, 50
|
||||
|
||||
H = np.eye(3)
|
||||
S = scale * np.sin(theta)
|
||||
C = scale * np.cos(theta)
|
||||
|
||||
H[:2, :2] = [[C, -S], [S, C]]
|
||||
H[:2, 2] = [tx, ty]
|
||||
|
||||
for mode in ('constant', 'mirror', 'wrap'):
|
||||
print 'Transform mode:', mode
|
||||
|
||||
p0 = homography(img, H, mode=mode, order=1)
|
||||
p1 = fast_homography(img, H, mode=mode)
|
||||
p1 = np.round(p1)
|
||||
|
||||
## import matplotlib.pyplot as plt
|
||||
## f, (ax0, ax1, ax2, ax3) = plt.subplots(1, 4)
|
||||
## ax0.imshow(img)
|
||||
## ax1.imshow(p0, cmap=plt.cm.gray)
|
||||
## ax2.imshow(p1, cmap=plt.cm.gray)
|
||||
## ax3.imshow(np.abs(p0 - p1), cmap=plt.cm.gray)
|
||||
## plt.show()
|
||||
|
||||
d = np.mean(np.abs(p0 - p1))
|
||||
assert d < 0.2
|
||||
|
||||
|
||||
if __name__ == "__main__":
|
||||
from numpy.testing import run_module_suite
|
||||
run_module_suite()
|
||||
|
||||
Reference in New Issue
Block a user