1 // This file is part of OpenMVG, an Open Multiple View Geometry C++ library. 2 3 // Copyright (c) 2014 Pierre MOULON. 4 5 // This Source Code Form is subject to the terms of the Mozilla Public 6 // License, v. 2.0. If a copy of the MPL was not distributed with this 7 // file, You can obtain one at http://mozilla.org/MPL/2.0/. 8 9 #ifndef OPENMVG_FEATURES_DIPOLE_DIPOLE_DESCRIPTOR_HPP 10 #define OPENMVG_FEATURES_DIPOLE_DIPOLE_DESCRIPTOR_HPP 11 12 13 #include "openMVG/features/feature.hpp" 14 #include "openMVG/features/descriptor.hpp" 15 #include "openMVG/image/image_container.hpp" 16 #include "openMVG/image/sample.hpp" 17 18 //------------------ 19 //-- Bibliography -- 20 //------------------ 21 //- [1] "New local descriptor based on dissociated dipoles" 22 //- Authors: Alexis Joly. 23 //- Date: December 2007. 24 //- Conference: CIVR. 25 26 namespace openMVG 27 { 28 namespace features 29 { 30 // Create a DIPOLE describer [1]. 31 // 32 // Note : 33 // - Angle is in radians. 34 // - data the output array (must be allocated to 20 values). 35 template<typename Real> PickNaiveDipole(const image::Image<Real> & image,float x,float y,float scale,float angle,float * data)36 void PickNaiveDipole 37 ( 38 const image::Image<Real> & image, 39 float x, 40 float y, 41 float scale, 42 float angle, 43 float * data 44 ) 45 { 46 // Use bilinear sampling 47 const image::Sampler2d<image::SamplerLinear> sampler; 48 // Setup the rotation center. 49 const float & cx = x, & cy = y; 50 51 const float lambda1 = scale; 52 const float lambda2 = lambda1 / 2.0f; 53 static const float angleSubdiv = 2.0f * M_PI / 12.0f; 54 const float memoizeCos[12] = 55 { 56 std::cos(angle + 0.f * angleSubdiv), 57 std::cos(angle + 1.f * angleSubdiv), 58 std::cos(angle + 2.f * angleSubdiv), 59 std::cos(angle + 3.f * angleSubdiv), 60 std::cos(angle + 4.f * angleSubdiv), 61 std::cos(angle + 5.f * angleSubdiv), 62 std::cos(angle + 6.f * angleSubdiv), 63 std::cos(angle + 7.f * angleSubdiv), 64 std::cos(angle + 8.f * angleSubdiv), 65 std::cos(angle + 9.f * angleSubdiv), 66 std::cos(angle + 10.f * angleSubdiv), 67 std::cos(angle + 11.f * angleSubdiv) 68 }; 69 const float memoizeSin[12] = 70 { 71 std::sin(angle + 0.f * angleSubdiv), 72 std::sin(angle + 1.f * angleSubdiv), 73 std::sin(angle + 2.f * angleSubdiv), 74 std::sin(angle + 3.f * angleSubdiv), 75 std::sin(angle + 4.f * angleSubdiv), 76 std::sin(angle + 5.f * angleSubdiv), 77 std::sin(angle + 6.f * angleSubdiv), 78 std::sin(angle + 7.f * angleSubdiv), 79 std::sin(angle + 8.f * angleSubdiv), 80 std::sin(angle + 9.f * angleSubdiv), 81 std::sin(angle + 10.f * angleSubdiv), 82 std::sin(angle + 11.f * angleSubdiv) 83 }; 84 85 Vecf dipoleF1(12); 86 for (int i = 0; i < 12; ++i) 87 { 88 const float xi = cx + lambda1 * memoizeCos[i]; 89 const float yi = cy + lambda1 * memoizeSin[i]; 90 dipoleF1(i) = sampler(image, yi, xi); 91 } 92 Matf A(8,12); 93 A << 0, 0, 0, 1, 0, 0, 0, 0, 0,-1, 0, 0, 94 0,-1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 95 0, 0, 0, 0, 0,-1, 0, 0, 0, 0, 0, 1, 96 0, 0, 0, 0, 1, 0, 0,-1, 0, 0, 0, 0, 97 0, 0, 0, 0, 0, 0, 1, 0, 0,-1, 0, 0, 98 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0,-1, 99 0,-1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 100 1, 0, 0,-1, 0, 0, 0, 0, 0, 0, 0, 0; 101 102 // Add the second order F2 dipole 103 Vecf dipoleF2(12); 104 for (int i = 0; i < 12; ++i) 105 { 106 const float xi = cx + (lambda1 + lambda2) * memoizeCos[i]; 107 const float yi = cy + (lambda1 + lambda2) * memoizeSin[i]; 108 109 const float xii = cx + (lambda1 - lambda2) * memoizeCos[i]; 110 const float yii = cy + (lambda1 - lambda2) * memoizeSin[i]; 111 112 // Bilinear interpolation 113 dipoleF2(i) = sampler(image, yi, xi) - sampler(image, yii, xii); 114 } 115 // Normalize to be affine luminance invariant (a*I(x,y)+b). 116 Map<Vecf> dataMap( data, 20); 117 dataMap.block<8,1>(0,0) = (A * dipoleF1).normalized(); 118 dataMap.block<12,1>(8,0) = dipoleF2.normalized(); 119 } 120 121 // Pick an angular smoothed dipole 122 template<typename Real> PickASDipole(const image::Image<Real> & image,float x,float y,float scale,float angle,float * data)123 void PickASDipole 124 ( 125 const image::Image<Real> & image, 126 float x, 127 float y, 128 float scale, 129 float angle, 130 float * data) 131 { 132 const image::Sampler2d<image::SamplerLinear> sampler; 133 // Setup the rotation center. 134 const float & cx = x, & cy = y; 135 136 const float lambda1 = scale; 137 const float lambda2 = lambda1 / 2.0f; 138 const float angleSubdiv = 2.0f * M_PI / 12.0f; 139 140 //-- First order dipole: 141 Vecf dipoleF1(12); 142 for (int i = 0; i < 12; ++i) 143 { 144 const float xi = cx + lambda1 * std::cos(angle + i * angleSubdiv); 145 const float yi = cy + lambda1 * std::sin(angle + i * angleSubdiv); 146 const float xi0 = cx + lambda1 * std::cos(angle + i * angleSubdiv - angleSubdiv/2.0); 147 const float yi0 = cy + lambda1 * std::sin(angle + i * angleSubdiv - angleSubdiv/2.0); 148 const float xi3 = cx + lambda1 * std::cos(angle + i * angleSubdiv + angleSubdiv/2.0); 149 const float yi3 = cy + lambda1 * std::sin(angle + i * angleSubdiv + angleSubdiv/2.0); 150 // Bilinear interpolation 151 dipoleF1(i) = 152 (sampler(image, yi, xi) + 153 sampler(image, yi0, xi0) + 154 sampler(image, yi3, xi3))/3.0f; 155 } 156 Matf A(8,12); 157 A << 0, 0, 0, 1, 0, 0, 0, 0, 0,-1, 0, 0, 158 0,-1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 159 0, 0, 0, 0, 0,-1, 0, 0, 0, 0, 0, 1, 160 0, 0, 0, 0, 1, 0, 0,-1, 0, 0, 0, 0, 161 0, 0, 0, 0, 0, 0, 1, 0, 0,-1, 0, 0, 162 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0,-1, 163 0,-1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 164 1, 0, 0,-1, 0, 0, 0, 0, 0, 0, 0, 0; 165 166 //-- Second order dipole: 167 Vecf dipoleF2(12); 168 for (int i = 0; i < 12; ++i) 169 { 170 const double angleSample = i * angleSubdiv; 171 const float xi = cx + (lambda1 + lambda2) * std::cos(angle + angleSample); 172 const float yi = cy + (lambda1 + lambda2) * std::sin(angle + angleSample); 173 const float xi0 = cx + (lambda1 + lambda2) * std::cos(angle + angleSample - angleSubdiv/2.0); 174 const float yi0 = cy + (lambda1 + lambda2) * std::sin(angle + angleSample - angleSubdiv/2.0); 175 const float xi3 = cx + (lambda1 + lambda2) * std::cos(angle + angleSample + angleSubdiv/2.0); 176 const float yi3 = cy + (lambda1 + lambda2) * std::sin(angle + angleSample + angleSubdiv/2.0); 177 178 const float xii = cx + (lambda1 - lambda2) * std::cos(angle + angleSample); 179 const float yii = cy + (lambda1 - lambda2) * std::sin(angle + angleSample); 180 const float xii0 = cx + (lambda1 - lambda2) * std::cos(angle + angleSample - angleSubdiv/2.0); 181 const float yii0 = cy + (lambda1 - lambda2) * std::sin(angle + angleSample - angleSubdiv/2.0); 182 const float xii3 = cx + (lambda1 - lambda2) * std::cos(angle + angleSample + angleSubdiv/2.0); 183 const float yii3 = cy + (lambda1 - lambda2) * std::sin(angle + angleSample + angleSubdiv/2.0); 184 185 // Bilinear interpolation 186 dipoleF2(i) = 187 (sampler(image, yi, xi) + 188 sampler(image, yi0, xi0) + 189 sampler(image, yi3, xi3)) /3.0f 190 - 191 (sampler(image, yii, xii) + 192 sampler(image, yii0, xii0) + 193 sampler(image, yii3, xii3)) /3.0f; 194 } 195 // Normalize to be affine luminance invariant (a*I(x,y)+b). 196 Map<Vecf> dataMap( data, 20); 197 dataMap.block<8,1>(0,0) = (A * dipoleF1).normalized(); 198 dataMap.block<12,1>(8,0) = dipoleF2.normalized(); 199 } 200 201 /** 202 ** @brief Compute DIPOLE descriptor for a given interest point 203 ** @param Li Input image 204 ** @param ipt Input interest point 205 ** @param desc output descriptor (floating point descriptor) 206 ** @param bAngularSmoothedDipole Tell if we must extract an upright or an 207 ** angular smoothed dipole 208 ** @param magnif_factor Scaling factor used to rescale the dipole sampling 209 **/ 210 template<typename Real> ComputeDipoleDescriptor(const image::Image<Real> & Li,const SIOPointFeature & ipt,Descriptor<float,20> & desc,bool bAngularSmoothedDipole=true,const float magnif_factor=3.5f)211 void ComputeDipoleDescriptor 212 ( 213 const image::Image<Real> & Li, 214 const SIOPointFeature & ipt, 215 Descriptor<float, 20> & desc, 216 bool bAngularSmoothedDipole = true, 217 const float magnif_factor = 3.5f 218 ) 219 { 220 if (bAngularSmoothedDipole) 221 PickASDipole(Li, ipt.x(), ipt.y(), ipt.scale() * magnif_factor, ipt.orientation(), &desc[0]); 222 else 223 PickNaiveDipole(Li, ipt.x(), ipt.y(), ipt.scale() * magnif_factor, ipt.orientation(), &desc[0]); 224 } 225 226 } // namespace features 227 } // namespace openMVG 228 229 #endif // OPENMVG_FEATURES_DIPOLE_DIPOLE_DESCRIPTOR_HPP 230