Line | Count | Source (jump to first uncovered line) |
1 | | /****************************************************************************** |
2 | | * |
3 | | * Project: MapServer |
4 | | * Purpose: Inverse Distance Weighted layer interpolation. |
5 | | * Author: Hermes L. Herrera Martinez and the MapServer team. |
6 | | * Thanks: Thomas Bonfort |
7 | | * |
8 | | ****************************************************************************** |
9 | | * Copyright (c) 2014 Regents of the University of Minnesota. |
10 | | * |
11 | | * Permission is hereby granted, free of charge, to any person obtaining a |
12 | | * copy of this software and associated documentation files (the "Software"), |
13 | | * to deal in the Software without restriction, including without limitation |
14 | | * the rights to use, copy, modify, merge, publish, distribute, sublicense, |
15 | | * and/or sell copies of the Software, and to permit persons to whom the |
16 | | * Software is furnished to do so, subject to the following conditions: |
17 | | * |
18 | | * The above copyright notice and this permission notice shall be included in |
19 | | * all copies of this Software or works derived from this Software. |
20 | | * |
21 | | * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS |
22 | | * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, |
23 | | * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL |
24 | | * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER |
25 | | * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING |
26 | | * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER |
27 | | * DEALINGS IN THE SOFTWARE. |
28 | | *****************************************************************************/ |
29 | | |
30 | | #include "mapserver.h" |
31 | | #include <float.h> |
32 | 0 | #define EPSILON 0.000000001 |
33 | | #include <time.h> |
34 | | |
35 | | void msIdw(float *xyz, int width, int height, int npoints, |
36 | | interpolationProcessingParams *interpParams, |
37 | 0 | unsigned char *iValues) { |
38 | 0 | int i, j, index; |
39 | 0 | int radius = interpParams->radius; |
40 | 0 | float power = interpParams->power; |
41 | 0 | for (j = 0; j < height; j++) { |
42 | 0 | for (i = 0; i < width; i++) { |
43 | 0 | double den = EPSILON, num = 0; |
44 | 0 | for (index = 0; index < npoints * 3; index += 3) { |
45 | 0 | double d = (xyz[index] - i) * (xyz[index] - i) + |
46 | 0 | (xyz[index + 1] - j) * (xyz[index + 1] - j); |
47 | 0 | if (radius * radius > d) { |
48 | 0 | double w = 1.0 / (pow(d, power) + EPSILON); |
49 | 0 | num += w * xyz[index + 2]; |
50 | 0 | den += w; |
51 | 0 | } |
52 | 0 | } |
53 | 0 | iValues[j * width + i] = num / den; |
54 | 0 | } |
55 | 0 | } |
56 | 0 | } |
57 | | |
58 | | void msIdwProcessing(layerObj *layer, |
59 | 0 | interpolationProcessingParams *interpParams) { |
60 | 0 | const char *interpParamsProcessing = |
61 | 0 | msLayerGetProcessingKey(layer, "IDW_POWER"); |
62 | 0 | if (interpParamsProcessing) { |
63 | 0 | interpParams->power = atof(interpParamsProcessing); |
64 | 0 | } else { |
65 | 0 | interpParams->power = 1.0; |
66 | 0 | } |
67 | |
|
68 | 0 | interpParamsProcessing = msLayerGetProcessingKey(layer, "IDW_RADIUS"); |
69 | 0 | if (interpParamsProcessing) { |
70 | 0 | interpParams->radius = atof(interpParamsProcessing); |
71 | 0 | } else { |
72 | 0 | interpParams->radius = MAX(layer->map->width, layer->map->height); |
73 | 0 | } |
74 | |
|
75 | 0 | interpParamsProcessing = |
76 | 0 | msLayerGetProcessingKey(layer, "IDW_COMPUTE_BORDERS"); |
77 | 0 | if (interpParamsProcessing && strcasecmp(interpParamsProcessing, "OFF")) { |
78 | 0 | interpParams->expand_searchrect = 1; |
79 | 0 | } else { |
80 | 0 | interpParams->expand_searchrect = 0; |
81 | 0 | } |
82 | 0 | } |