001/** 002 * Copyright (c) 2011, The University of Southampton and the individual contributors. 003 * All rights reserved. 004 * 005 * Redistribution and use in source and binary forms, with or without modification, 006 * are permitted provided that the following conditions are met: 007 * 008 * * Redistributions of source code must retain the above copyright notice, 009 * this list of conditions and the following disclaimer. 010 * 011 * * Redistributions in binary form must reproduce the above copyright notice, 012 * this list of conditions and the following disclaimer in the documentation 013 * and/or other materials provided with the distribution. 014 * 015 * * Neither the name of the University of Southampton nor the names of its 016 * contributors may be used to endorse or promote products derived from this 017 * software without specific prior written permission. 018 * 019 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND 020 * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED 021 * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE 022 * DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR 023 * ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES 024 * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; 025 * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON 026 * ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT 027 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS 028 * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 029 */ 030package org.openimaj.image.feature.local.descriptor.gradient; 031 032import org.openimaj.citation.annotation.Reference; 033import org.openimaj.citation.annotation.ReferenceType; 034import org.openimaj.citation.annotation.References; 035import org.openimaj.feature.OrientedFeatureVector; 036import org.openimaj.util.array.ArrayUtils; 037 038/** 039 * An extractor for SIFT features. SIFT features are basically multiple edge 040 * orientation histograms constructed over a spatial grid. Samples added to the 041 * SIFT histogram are weighted with a Gaussian based on the distance from the 042 * centre of the sampling window. Samples are also blurred across histogram bins 043 * in both the spatial and orientation directions. 044 * 045 * Based on Section 6 of Lowe's IJCV paper 046 * 047 * @author Jonathon Hare (jsh2@ecs.soton.ac.uk) 048 */ 049@References(references = { 050 @Reference( 051 type = ReferenceType.Article, 052 author = { "David Lowe" }, 053 title = "Distinctive image features from scale-invariant keypoints", 054 year = "2004", 055 journal = "IJCV", 056 pages = { "91", "110" }, 057 month = "January", 058 number = "2", 059 volume = "60"), 060 @Reference( 061 type = ReferenceType.Inproceedings, 062 author = { "David Lowe" }, 063 title = "Object recognition from local scale-invariant features", 064 year = "1999", 065 booktitle = "Proc. of the International Conference on Computer Vision {ICCV}", 066 pages = { "1150", "1157" } 067 ) 068}) 069public class SIFTFeatureProvider implements GradientFeatureProvider, GradientFeatureProviderFactory { 070 private final static float TWO_PI_FLOAT = (float) (Math.PI * 2); 071 072 /** Number of orientation bins in the histograms */ 073 protected int numOriBins = 8; 074 075 /** Number of spatial bins for the x and y directions in the histograms */ 076 protected int numSpatialBins = 4; 077 078 /** Threshold for the maximum allowed value in the histogram */ 079 protected float valueThreshold = 0.2f; 080 081 /** 082 * 2 times the weighting Gaussian squared (normalised to the patch size in 083 * terms of spatial bins) 084 */ 085 protected float sigmaSq2 = 0.5f; 086 087 protected float gaussianSigma = 1; 088 089 protected float[] vec; 090 091 protected float patchOrientation; 092 093 /** 094 * Construct a {@link SIFTFeatureProvider} with the default parameters. 095 */ 096 public SIFTFeatureProvider() { 097 this.vec = new float[numSpatialBins * numSpatialBins * numOriBins]; 098 } 099 100 /** 101 * Construct a {@link SIFTFeatureProvider} with the provided options. 102 * 103 * @param numOriBins 104 * the number of orientation bins (default 8) 105 * @param numSpatialBins 106 * the number of spatial bins in each direction (default 4) 107 */ 108 public SIFTFeatureProvider(int numOriBins, int numSpatialBins) { 109 this(numOriBins, numSpatialBins, 0.2f, 1.0f); 110 } 111 112 /** 113 * Construct a {@link SIFTFeatureProvider} with the provided options. 114 * 115 * @param numOriBins 116 * the number of orientation bins (default 8) 117 * @param numSpatialBins 118 * the number of spatial bins in each direction (default 4) 119 * @param valueThreshold 120 * threshold for the maximum value allowed in the histogram 121 * (default 0.2) 122 * @param gaussianSigma 123 * the width of the Gaussian used for weighting samples, relative 124 * to the half-width of the sampling window (default 1.0). 125 */ 126 public SIFTFeatureProvider(int numOriBins, int numSpatialBins, float valueThreshold, float gaussianSigma) { 127 this.numOriBins = numOriBins; 128 this.numSpatialBins = numSpatialBins; 129 this.valueThreshold = valueThreshold; 130 this.gaussianSigma = gaussianSigma; 131 this.vec = new float[numSpatialBins * numSpatialBins * numOriBins]; 132 133 // calculate the variance of the Gaussian weighting 134 final float sigma = gaussianSigma / (0.5f * numSpatialBins); // indexSigma 135 // is 136 // proportional 137 // to 138 // half 139 // the 140 // index 141 // size 142 sigmaSq2 = 2 * sigma * sigma; 143 } 144 145 @Override 146 public void addSample(float x, float y, float gradmag, float gradori) { 147 // calculate weight based on Gaussian at the centre: 148 final float dx = 0.5f - x; 149 final float dy = 0.5f - y; 150 final float weight = (float) Math.exp(-(dx * dx + dy * dy) / sigmaSq2); 151 152 // weight the magnitude 153 final float wmag = weight * gradmag; 154 155 // adjust the gradient angle to be relative to the patch angle 156 float ori = gradori - patchOrientation; 157 158 // adjust range to 0<=ori<2PI 159 ori = ((ori %= TWO_PI_FLOAT) >= 0 ? ori : (ori + TWO_PI_FLOAT)); 160 161 // now add the sample to the correct bins 162 interpolateSample(x, y, wmag, ori); 163 } 164 165 /** 166 * Spread the sample around the closest bins in the histogram. If there are 167 * 4 spatial bins in each direction, then a sample at 0.25 would get added 168 * equally to bins [0] and [1] in the x-direction. 169 * 170 * @param x 171 * the normalised x-coordinate 172 * @param y 173 * the normalised y-coordinate 174 * @param magnitude 175 * the magnitude of the sample 176 * @param orientation 177 * the angle of the sample 178 */ 179 protected void interpolateSample(float x, float y, float magnitude, float orientation) { 180 final float px = numSpatialBins * x - 0.5f; // px is now 181 // 0.5<=px<=indexSize-0.5 182 final float py = numSpatialBins * y - 0.5f; // py is now 183 // 0.5<=py<=indexSize-0.5 184 final float po = numOriBins * orientation / TWO_PI_FLOAT; // po is now 185 // 0<=po<oriSize 186 187 // the integer parts - corresponding to the left (or equivalent) bin 188 // that the sample falls in 189 final int xi = (int) Math.floor(px); 190 final int yi = (int) Math.floor(py); 191 final int oi = (int) Math.floor(po); 192 193 // the fractional parts - corresponding to how much goes in the right 194 // bin 195 // 1-xf goes in the left bin 196 final float xf = px - xi; 197 final float yf = py - yi; 198 final float of = po - oi; 199 200 // now spread the sample around a 2x2x2 cube (left bin, right bin each 201 // each dim 202 // + combinations) 203 for (int yy = 0; yy < 2; yy++) { 204 final int yindex = yi + yy; 205 206 if (yindex >= 0 && yindex < numSpatialBins) { 207 final float yweight = magnitude * ((yy == 0) ? 1.0f - yf : yf); 208 209 for (int xx = 0; xx < 2; xx++) { 210 final int xindex = xi + xx; 211 212 if (xindex >= 0 && xindex < numSpatialBins) { 213 final float xweight = yweight * ((xx == 0) ? 1.0f - xf : xf); 214 215 for (int oo = 0; oo < 2; oo++) { 216 int oindex = oi + oo; 217 218 if (oindex >= numOriBins) 219 oindex = 0; // Orientation wraps at 2PI. 220 221 final float oweight = xweight * ((oo == 0) ? 1.0f - of : of); 222 223 vec[(numSpatialBins * numOriBins * yindex) + (numOriBins * xindex) + oindex] += oweight; 224 } 225 } 226 } 227 } 228 } 229 } 230 231 @Override 232 public OrientedFeatureVector getFeatureVector() { 233 ArrayUtils.normalise(vec); 234 235 boolean changed = false; 236 for (int i = 0; i < vec.length; i++) { 237 if (vec[i] > valueThreshold) { 238 vec[i] = valueThreshold; 239 changed = true; 240 } 241 } 242 243 if (changed) 244 ArrayUtils.normalise(vec); 245 246 // Construct the actual feature vector 247 final OrientedFeatureVector fv = new OrientedFeatureVector(vec.length, patchOrientation); 248 for (int i = 0; i < vec.length; i++) { 249 final int intval = (int) (512.0 * vec[i]); 250 251 fv.values[i] = (byte) (Math.min(255, intval) - 128); 252 } 253 254 return fv; 255 } 256 257 @Override 258 public void setPatchOrientation(float patchOrientation) { 259 this.patchOrientation = patchOrientation; 260 } 261 262 @Override 263 public GradientFeatureProvider newProvider() { 264 return new SIFTFeatureProvider(numOriBins, numSpatialBins, valueThreshold, gaussianSigma); 265 } 266 267 @Override 268 public float getOversamplingAmount() { 269 // in order to ensure a smooth interpolation, 270 // half a bin's width needs to be added to the 271 // sampling region all the way around the 272 // sampling square (so edge bins get partial 273 // contributions from outside the square) 274 return 1.0f / numSpatialBins / 2; 275 } 276}