1 /*
2 * Licensed to the Apache Software Foundation (ASF) under one or more
3 * contributor license agreements. See the NOTICE file distributed with
4 * this work for additional information regarding copyright ownership.
5 * The ASF licenses this file to You under the Apache License, Version 2.0
6 * (the "License"); you may not use this file except in compliance with
7 * the License. You may obtain a copy of the License at
8 *
9 * https://www.apache.org/licenses/LICENSE-2.0
10 *
11 * Unless required by applicable law or agreed to in writing, software
12 * distributed under the License is distributed on an "AS IS" BASIS,
13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 * See the License for the specific language governing permissions and
15 * limitations under the License.
16 */
17
18 /*
19 * This is not the original file distributed by the Apache Software Foundation
20 * It has been modified by the Hipparchus project
21 */
22 package org.hipparchus.analysis.interpolation;
23
24 import java.util.ArrayList;
25 import java.util.List;
26
27 import org.hipparchus.exception.LocalizedCoreFormats;
28 import org.hipparchus.exception.MathIllegalArgumentException;
29 import org.hipparchus.exception.MathIllegalStateException;
30 import org.hipparchus.random.UnitSphereRandomVectorGenerator;
31 import org.hipparchus.util.FastMath;
32 import org.hipparchus.util.MathArrays;
33 import org.hipparchus.util.MathUtils;
34
35 /**
36 * Utility class for the {@link MicrosphereProjectionInterpolator} algorithm.
37 *
38 */
39 public class InterpolatingMicrosphere {
40 /** Microsphere. */
41 private final List<Facet> microsphere;
42 /** Microsphere data. */
43 private final List<FacetData> microsphereData;
44 /** Space dimension. */
45 private final int dimension;
46 /** Number of surface elements. */
47 private final int size;
48 /** Maximum fraction of the facets that can be dark. */
49 private final double maxDarkFraction;
50 /** Lowest non-zero illumination. */
51 private final double darkThreshold;
52 /** Background value. */
53 private final double background;
54
55 /**
56 * Create an unitialiazed sphere.
57 * Sub-classes are responsible for calling the {@code add(double[]) add}
58 * method in order to initialize all the sphere's facets.
59 *
60 * @param dimension Dimension of the data space.
61 * @param size Number of surface elements of the sphere.
62 * @param maxDarkFraction Maximum fraction of the facets that can be dark.
63 * If the fraction of "non-illuminated" facets is larger, no estimation
64 * of the value will be performed, and the {@code background} value will
65 * be returned instead.
66 * @param darkThreshold Value of the illumination below which a facet is
67 * considered dark.
68 * @param background Value returned when the {@code maxDarkFraction}
69 * threshold is exceeded.
70 * @throws MathIllegalArgumentException if {@code dimension <= 0}
71 * or {@code size <= 0}.
72 * @throws MathIllegalArgumentException if {@code darkThreshold < 0}.
73 * @throws MathIllegalArgumentException if {@code maxDarkFraction} does not
74 * belong to the interval {@code [0, 1]}.
75 */
76 protected InterpolatingMicrosphere(int dimension,
77 int size,
78 double maxDarkFraction,
79 double darkThreshold,
80 double background) {
81 if (dimension <= 0) {
82 throw new MathIllegalArgumentException(LocalizedCoreFormats.NUMBER_TOO_SMALL_BOUND_EXCLUDED,
83 dimension, 0);
84 }
85 if (size <= 0) {
86 throw new MathIllegalArgumentException(LocalizedCoreFormats.NUMBER_TOO_SMALL_BOUND_EXCLUDED,
87 size, 0);
88 }
89 MathUtils.checkRangeInclusive(maxDarkFraction, 0, 1);
90 if (darkThreshold < 0) {
91 throw new MathIllegalArgumentException(LocalizedCoreFormats.NUMBER_TOO_SMALL, darkThreshold, 0);
92 }
93
94 this.dimension = dimension;
95 this.size = size;
96 this.maxDarkFraction = maxDarkFraction;
97 this.darkThreshold = darkThreshold;
98 this.background = background;
99 microsphere = new ArrayList<>(size);
100 microsphereData = new ArrayList<>(size);
101 }
102
103 /**
104 * Create a sphere from randomly sampled vectors.
105 *
106 * @param dimension Dimension of the data space.
107 * @param size Number of surface elements of the sphere.
108 * @param rand Unit vector generator for creating the microsphere.
109 * @param maxDarkFraction Maximum fraction of the facets that can be dark.
110 * If the fraction of "non-illuminated" facets is larger, no estimation
111 * of the value will be performed, and the {@code background} value will
112 * be returned instead.
113 * @param darkThreshold Value of the illumination below which a facet
114 * is considered dark.
115 * @param background Value returned when the {@code maxDarkFraction}
116 * threshold is exceeded.
117 * @throws MathIllegalArgumentException if the size of the generated
118 * vectors does not match the dimension set in the constructor.
119 * @throws MathIllegalArgumentException if {@code dimension <= 0}
120 * or {@code size <= 0}.
121 * @throws MathIllegalArgumentException if {@code darkThreshold < 0}.
122 * @throws MathIllegalArgumentException if {@code maxDarkFraction} does not
123 * belong to the interval {@code [0, 1]}.
124 */
125 public InterpolatingMicrosphere(int dimension,
126 int size,
127 double maxDarkFraction,
128 double darkThreshold,
129 double background,
130 UnitSphereRandomVectorGenerator rand) {
131 this(dimension, size, maxDarkFraction, darkThreshold, background);
132
133 // Generate the microsphere normals, assuming that a number of
134 // randomly generated normals will represent a sphere.
135 for (int i = 0; i < size; i++) {
136 add(rand.nextVector(), false);
137 }
138 }
139
140 /**
141 * Copy constructor.
142 *
143 * @param other Instance to copy.
144 */
145 protected InterpolatingMicrosphere(InterpolatingMicrosphere other) {
146 dimension = other.dimension;
147 size = other.size;
148 maxDarkFraction = other.maxDarkFraction;
149 darkThreshold = other.darkThreshold;
150 background = other.background;
151
152 // Field can be shared.
153 microsphere = other.microsphere;
154
155 // Field must be copied.
156 microsphereData = new ArrayList<>(size);
157 for (FacetData fd : other.microsphereData) {
158 microsphereData.add(new FacetData(fd.illumination(), fd.sample()));
159 }
160 }
161
162 /**
163 * Perform a copy.
164 *
165 * @return a copy of this instance.
166 */
167 public InterpolatingMicrosphere copy() {
168 return new InterpolatingMicrosphere(this);
169 }
170
171 /**
172 * Get the space dimensionality.
173 *
174 * @return the number of space dimensions.
175 */
176 public int getDimension() {
177 return dimension;
178 }
179
180 /**
181 * Get the size of the sphere.
182 *
183 * @return the number of surface elements of the microspshere.
184 */
185 public int getSize() {
186 return size;
187 }
188
189 /**
190 * Estimate the value at the requested location.
191 * This microsphere is placed at the given {@code point}, contribution
192 * of the given {@code samplePoints} to each sphere facet is computed
193 * (illumination) and the interpolation is performed (integration of
194 * the illumination).
195 *
196 * @param point Interpolation point.
197 * @param samplePoints Sampling data points.
198 * @param sampleValues Sampling data values at the corresponding
199 * {@code samplePoints}.
200 * @param exponent Exponent used in the power law that computes
201 * the weights (distance dimming factor) of the sample data.
202 * @param noInterpolationTolerance When the distance between the
203 * {@code point} and one of the {@code samplePoints} is less than
204 * this value, no interpolation will be performed, and the value
205 * of the sample will just be returned.
206 * @return the estimated value at the given {@code point}.
207 * @throws MathIllegalArgumentException if {@code exponent < 0}.
208 */
209 public double value(double[] point,
210 double[][] samplePoints,
211 double[] sampleValues,
212 double exponent,
213 double noInterpolationTolerance) {
214 if (exponent < 0) {
215 throw new MathIllegalArgumentException(LocalizedCoreFormats.NUMBER_TOO_SMALL, exponent, 0);
216 }
217
218 clear();
219
220 // Contribution of each sample point to the illumination of the
221 // microsphere's facets.
222 final int numSamples = samplePoints.length;
223 for (int i = 0; i < numSamples; i++) {
224 // Vector between interpolation point and current sample point.
225 final double[] diff = MathArrays.ebeSubtract(samplePoints[i], point);
226 final double diffNorm = MathArrays.safeNorm(diff);
227
228 if (FastMath.abs(diffNorm) < noInterpolationTolerance) {
229 // No need to interpolate, as the interpolation point is
230 // actually (very close to) one of the sampled points.
231 return sampleValues[i];
232 }
233
234 final double weight = FastMath.pow(diffNorm, -exponent);
235 illuminate(diff, sampleValues[i], weight);
236 }
237
238 return interpolate();
239 }
240
241 /**
242 * Replace {@code i}-th facet of the microsphere.
243 * Method for initializing the microsphere facets.
244 *
245 * @param normal Facet's normal vector.
246 * @param copy Whether to copy the given array.
247 * @throws MathIllegalArgumentException if the length of {@code n}
248 * does not match the space dimension.
249 * @throws MathIllegalStateException if the method has been called
250 * more times than the size of the sphere.
251 */
252 protected void add(double[] normal,
253 boolean copy) {
254 if (microsphere.size() >= size) {
255 throw new MathIllegalStateException(LocalizedCoreFormats.MAX_COUNT_EXCEEDED, size);
256 }
257 if (normal.length > dimension) {
258 throw new MathIllegalArgumentException(LocalizedCoreFormats.DIMENSIONS_MISMATCH,
259 normal.length, dimension);
260 }
261
262 microsphere.add(new Facet(copy ? normal.clone() : normal));
263 microsphereData.add(new FacetData(0d, 0d));
264 }
265
266 /**
267 * Interpolation.
268 *
269 * @return the value estimated from the current illumination of the
270 * microsphere.
271 */
272 private double interpolate() {
273 // Number of non-illuminated facets.
274 int darkCount = 0;
275
276 double value = 0;
277 double totalWeight = 0;
278 for (FacetData fd : microsphereData) {
279 final double iV = fd.illumination();
280 if (iV != 0d) {
281 value += iV * fd.sample();
282 totalWeight += iV;
283 } else {
284 ++darkCount;
285 }
286 }
287
288 final double darkFraction = darkCount / (double) size;
289
290 return darkFraction <= maxDarkFraction ?
291 value / totalWeight :
292 background;
293 }
294
295 /**
296 * Illumination.
297 *
298 * @param sampleDirection Vector whose origin is at the interpolation
299 * point and tail is at the sample location.
300 * @param sampleValue Data value of the sample.
301 * @param weight Weight.
302 */
303 private void illuminate(double[] sampleDirection,
304 double sampleValue,
305 double weight) {
306 for (int i = 0; i < size; i++) {
307 final double[] n = microsphere.get(i).getNormal();
308 final double cos = MathArrays.cosAngle(n, sampleDirection);
309
310 if (cos > 0) {
311 final double illumination = cos * weight;
312
313 if (illumination > darkThreshold &&
314 illumination > microsphereData.get(i).illumination()) {
315 microsphereData.set(i, new FacetData(illumination, sampleValue));
316 }
317 }
318 }
319 }
320
321 /**
322 * Reset the all the {@link Facet facets} data to zero.
323 */
324 private void clear() {
325 for (int i = 0; i < size; i++) {
326 microsphereData.set(i, new FacetData(0d, 0d));
327 }
328 }
329
330 /**
331 * Microsphere "facet" (surface element).
332 */
333 private static class Facet {
334 /** Normal vector characterizing a surface element. */
335 private final double[] normal;
336
337 /**
338 * @param n Normal vector characterizing a surface element
339 * of the microsphere. No copy is made.
340 */
341 Facet(double[] n) {
342 normal = n; // NOPMD - array cloning is taken care of at call site
343 }
344
345 /**
346 * Return a reference to the vector normal to this facet.
347 *
348 * @return the normal vector.
349 */
350 public double[] getNormal() {
351 return normal; // NOPMD - returning an internal array is intentional and documented here
352 }
353 }
354
355 /**
356 * Data associated with each {@link Facet}.
357 */
358 private static class FacetData {
359 /** Illumination received from the sample. */
360 private final double illumination;
361 /** Data value of the sample. */
362 private final double sample;
363
364 /**
365 * @param illumination Illumination.
366 * @param sample Data value.
367 */
368 FacetData(double illumination, double sample) {
369 this.illumination = illumination;
370 this.sample = sample;
371 }
372
373 /**
374 * Get the illumination.
375 * @return the illumination.
376 */
377 public double illumination() {
378 return illumination;
379 }
380
381 /**
382 * Get the data value.
383 * @return the data value.
384 */
385 public double sample() {
386 return sample;
387 }
388 }
389 }