comparison libs/commons-math-2.1/docs/apidocs/src-html/org/apache/commons/math/analysis/interpolation/SmoothingBicubicSplineInterpolator.html @ 13:cbf34dd4d7e6

commons-math-2.1 added
author dwinter
date Tue, 04 Jan 2011 10:02:07 +0100
parents
children
comparison
equal deleted inserted replaced
12:970d26a94fb7 13:cbf34dd4d7e6
1 <HTML>
2 <BODY BGCOLOR="white">
3 <PRE>
4 <FONT color="green">001</FONT> /*<a name="line.1"></a>
5 <FONT color="green">002</FONT> * Licensed to the Apache Software Foundation (ASF) under one or more<a name="line.2"></a>
6 <FONT color="green">003</FONT> * contributor license agreements. See the NOTICE file distributed with<a name="line.3"></a>
7 <FONT color="green">004</FONT> * this work for additional information regarding copyright ownership.<a name="line.4"></a>
8 <FONT color="green">005</FONT> * The ASF licenses this file to You under the Apache License, Version 2.0<a name="line.5"></a>
9 <FONT color="green">006</FONT> * (the "License"); you may not use this file except in compliance with<a name="line.6"></a>
10 <FONT color="green">007</FONT> * the License. You may obtain a copy of the License at<a name="line.7"></a>
11 <FONT color="green">008</FONT> *<a name="line.8"></a>
12 <FONT color="green">009</FONT> * http://www.apache.org/licenses/LICENSE-2.0<a name="line.9"></a>
13 <FONT color="green">010</FONT> *<a name="line.10"></a>
14 <FONT color="green">011</FONT> * Unless required by applicable law or agreed to in writing, software<a name="line.11"></a>
15 <FONT color="green">012</FONT> * distributed under the License is distributed on an "AS IS" BASIS,<a name="line.12"></a>
16 <FONT color="green">013</FONT> * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.<a name="line.13"></a>
17 <FONT color="green">014</FONT> * See the License for the specific language governing permissions and<a name="line.14"></a>
18 <FONT color="green">015</FONT> * limitations under the License.<a name="line.15"></a>
19 <FONT color="green">016</FONT> */<a name="line.16"></a>
20 <FONT color="green">017</FONT> package org.apache.commons.math.analysis.interpolation;<a name="line.17"></a>
21 <FONT color="green">018</FONT> <a name="line.18"></a>
22 <FONT color="green">019</FONT> import org.apache.commons.math.DimensionMismatchException;<a name="line.19"></a>
23 <FONT color="green">020</FONT> import org.apache.commons.math.MathRuntimeException;<a name="line.20"></a>
24 <FONT color="green">021</FONT> import org.apache.commons.math.MathException;<a name="line.21"></a>
25 <FONT color="green">022</FONT> import org.apache.commons.math.util.MathUtils;<a name="line.22"></a>
26 <FONT color="green">023</FONT> import org.apache.commons.math.analysis.UnivariateRealFunction;<a name="line.23"></a>
27 <FONT color="green">024</FONT> import org.apache.commons.math.analysis.BivariateRealFunction;<a name="line.24"></a>
28 <FONT color="green">025</FONT> import org.apache.commons.math.analysis.polynomials.PolynomialSplineFunction;<a name="line.25"></a>
29 <FONT color="green">026</FONT> <a name="line.26"></a>
30 <FONT color="green">027</FONT> /**<a name="line.27"></a>
31 <FONT color="green">028</FONT> * Generates a bicubic interpolation function.<a name="line.28"></a>
32 <FONT color="green">029</FONT> * Before interpolating, smoothing of the input data is performed using<a name="line.29"></a>
33 <FONT color="green">030</FONT> * splines.<a name="line.30"></a>
34 <FONT color="green">031</FONT> * See &lt;b&gt;Handbook on splines for the user&lt;/b&gt;, ISBN 084939404X,<a name="line.31"></a>
35 <FONT color="green">032</FONT> * chapter 2.<a name="line.32"></a>
36 <FONT color="green">033</FONT> *<a name="line.33"></a>
37 <FONT color="green">034</FONT> * @version $Revision$ $Date$<a name="line.34"></a>
38 <FONT color="green">035</FONT> * @since 2.1<a name="line.35"></a>
39 <FONT color="green">036</FONT> */<a name="line.36"></a>
40 <FONT color="green">037</FONT> public class SmoothingBicubicSplineInterpolator<a name="line.37"></a>
41 <FONT color="green">038</FONT> implements BivariateRealGridInterpolator {<a name="line.38"></a>
42 <FONT color="green">039</FONT> /**<a name="line.39"></a>
43 <FONT color="green">040</FONT> * {@inheritDoc}<a name="line.40"></a>
44 <FONT color="green">041</FONT> */<a name="line.41"></a>
45 <FONT color="green">042</FONT> public BivariateRealFunction interpolate(final double[] xval,<a name="line.42"></a>
46 <FONT color="green">043</FONT> final double[] yval,<a name="line.43"></a>
47 <FONT color="green">044</FONT> final double[][] zval)<a name="line.44"></a>
48 <FONT color="green">045</FONT> throws MathException, IllegalArgumentException {<a name="line.45"></a>
49 <FONT color="green">046</FONT> if (xval.length == 0 || yval.length == 0 || zval.length == 0) {<a name="line.46"></a>
50 <FONT color="green">047</FONT> throw MathRuntimeException.createIllegalArgumentException("no data");<a name="line.47"></a>
51 <FONT color="green">048</FONT> }<a name="line.48"></a>
52 <FONT color="green">049</FONT> if (xval.length != zval.length) {<a name="line.49"></a>
53 <FONT color="green">050</FONT> throw new DimensionMismatchException(xval.length, zval.length);<a name="line.50"></a>
54 <FONT color="green">051</FONT> }<a name="line.51"></a>
55 <FONT color="green">052</FONT> <a name="line.52"></a>
56 <FONT color="green">053</FONT> MathUtils.checkOrder(xval, 1, true);<a name="line.53"></a>
57 <FONT color="green">054</FONT> MathUtils.checkOrder(yval, 1, true);<a name="line.54"></a>
58 <FONT color="green">055</FONT> <a name="line.55"></a>
59 <FONT color="green">056</FONT> final int xLen = xval.length;<a name="line.56"></a>
60 <FONT color="green">057</FONT> final int yLen = yval.length;<a name="line.57"></a>
61 <FONT color="green">058</FONT> <a name="line.58"></a>
62 <FONT color="green">059</FONT> // Samples (first index is y-coordinate, i.e. subarray variable is x)<a name="line.59"></a>
63 <FONT color="green">060</FONT> // 0 &lt;= i &lt; xval.length<a name="line.60"></a>
64 <FONT color="green">061</FONT> // 0 &lt;= j &lt; yval.length<a name="line.61"></a>
65 <FONT color="green">062</FONT> // zX[j][i] = f(xval[i], yval[j])<a name="line.62"></a>
66 <FONT color="green">063</FONT> final double[][] zX = new double[yLen][xLen];<a name="line.63"></a>
67 <FONT color="green">064</FONT> for (int i = 0; i &lt; xLen; i++) {<a name="line.64"></a>
68 <FONT color="green">065</FONT> if (zval[i].length != yLen) {<a name="line.65"></a>
69 <FONT color="green">066</FONT> throw new DimensionMismatchException(zval[i].length, yLen);<a name="line.66"></a>
70 <FONT color="green">067</FONT> }<a name="line.67"></a>
71 <FONT color="green">068</FONT> <a name="line.68"></a>
72 <FONT color="green">069</FONT> for (int j = 0; j &lt; yLen; j++) {<a name="line.69"></a>
73 <FONT color="green">070</FONT> zX[j][i] = zval[i][j];<a name="line.70"></a>
74 <FONT color="green">071</FONT> }<a name="line.71"></a>
75 <FONT color="green">072</FONT> }<a name="line.72"></a>
76 <FONT color="green">073</FONT> <a name="line.73"></a>
77 <FONT color="green">074</FONT> final SplineInterpolator spInterpolator = new SplineInterpolator();<a name="line.74"></a>
78 <FONT color="green">075</FONT> <a name="line.75"></a>
79 <FONT color="green">076</FONT> // For each line y[j] (0 &lt;= j &lt; yLen), construct a 1D spline with<a name="line.76"></a>
80 <FONT color="green">077</FONT> // respect to variable x<a name="line.77"></a>
81 <FONT color="green">078</FONT> final PolynomialSplineFunction[] ySplineX = new PolynomialSplineFunction[yLen];<a name="line.78"></a>
82 <FONT color="green">079</FONT> for (int j = 0; j &lt; yLen; j++) {<a name="line.79"></a>
83 <FONT color="green">080</FONT> ySplineX[j] = spInterpolator.interpolate(xval, zX[j]);<a name="line.80"></a>
84 <FONT color="green">081</FONT> }<a name="line.81"></a>
85 <FONT color="green">082</FONT> <a name="line.82"></a>
86 <FONT color="green">083</FONT> // For every knot (xval[i], yval[j]) of the grid, calculate corrected<a name="line.83"></a>
87 <FONT color="green">084</FONT> // values zY_1<a name="line.84"></a>
88 <FONT color="green">085</FONT> final double[][] zY_1 = new double[xLen][yLen];<a name="line.85"></a>
89 <FONT color="green">086</FONT> for (int j = 0; j &lt; yLen; j++) {<a name="line.86"></a>
90 <FONT color="green">087</FONT> final PolynomialSplineFunction f = ySplineX[j];<a name="line.87"></a>
91 <FONT color="green">088</FONT> for (int i = 0; i &lt; xLen; i++) {<a name="line.88"></a>
92 <FONT color="green">089</FONT> zY_1[i][j] = f.value(xval[i]);<a name="line.89"></a>
93 <FONT color="green">090</FONT> }<a name="line.90"></a>
94 <FONT color="green">091</FONT> }<a name="line.91"></a>
95 <FONT color="green">092</FONT> <a name="line.92"></a>
96 <FONT color="green">093</FONT> // For each line x[i] (0 &lt;= i &lt; xLen), construct a 1D spline with<a name="line.93"></a>
97 <FONT color="green">094</FONT> // respect to variable y generated by array zY_1[i]<a name="line.94"></a>
98 <FONT color="green">095</FONT> final PolynomialSplineFunction[] xSplineY = new PolynomialSplineFunction[xLen];<a name="line.95"></a>
99 <FONT color="green">096</FONT> for (int i = 0; i &lt; xLen; i++) {<a name="line.96"></a>
100 <FONT color="green">097</FONT> xSplineY[i] = spInterpolator.interpolate(yval, zY_1[i]);<a name="line.97"></a>
101 <FONT color="green">098</FONT> }<a name="line.98"></a>
102 <FONT color="green">099</FONT> <a name="line.99"></a>
103 <FONT color="green">100</FONT> // For every knot (xval[i], yval[j]) of the grid, calculate corrected<a name="line.100"></a>
104 <FONT color="green">101</FONT> // values zY_2<a name="line.101"></a>
105 <FONT color="green">102</FONT> final double[][] zY_2 = new double[xLen][yLen];<a name="line.102"></a>
106 <FONT color="green">103</FONT> for (int i = 0; i &lt; xLen; i++) {<a name="line.103"></a>
107 <FONT color="green">104</FONT> final PolynomialSplineFunction f = xSplineY[i];<a name="line.104"></a>
108 <FONT color="green">105</FONT> for (int j = 0; j &lt; yLen; j++) {<a name="line.105"></a>
109 <FONT color="green">106</FONT> zY_2[i][j] = f.value(yval[j]);<a name="line.106"></a>
110 <FONT color="green">107</FONT> }<a name="line.107"></a>
111 <FONT color="green">108</FONT> }<a name="line.108"></a>
112 <FONT color="green">109</FONT> <a name="line.109"></a>
113 <FONT color="green">110</FONT> // Partial derivatives with respect to x at the grid knots<a name="line.110"></a>
114 <FONT color="green">111</FONT> final double[][] dZdX = new double[xLen][yLen];<a name="line.111"></a>
115 <FONT color="green">112</FONT> for (int j = 0; j &lt; yLen; j++) {<a name="line.112"></a>
116 <FONT color="green">113</FONT> final UnivariateRealFunction f = ySplineX[j].derivative();<a name="line.113"></a>
117 <FONT color="green">114</FONT> for (int i = 0; i &lt; xLen; i++) {<a name="line.114"></a>
118 <FONT color="green">115</FONT> dZdX[i][j] = f.value(xval[i]);<a name="line.115"></a>
119 <FONT color="green">116</FONT> }<a name="line.116"></a>
120 <FONT color="green">117</FONT> }<a name="line.117"></a>
121 <FONT color="green">118</FONT> <a name="line.118"></a>
122 <FONT color="green">119</FONT> // Partial derivatives with respect to y at the grid knots<a name="line.119"></a>
123 <FONT color="green">120</FONT> final double[][] dZdY = new double[xLen][yLen];<a name="line.120"></a>
124 <FONT color="green">121</FONT> for (int i = 0; i &lt; xLen; i++) {<a name="line.121"></a>
125 <FONT color="green">122</FONT> final UnivariateRealFunction f = xSplineY[i].derivative();<a name="line.122"></a>
126 <FONT color="green">123</FONT> for (int j = 0; j &lt; yLen; j++) {<a name="line.123"></a>
127 <FONT color="green">124</FONT> dZdY[i][j] = f.value(yval[j]);<a name="line.124"></a>
128 <FONT color="green">125</FONT> }<a name="line.125"></a>
129 <FONT color="green">126</FONT> }<a name="line.126"></a>
130 <FONT color="green">127</FONT> <a name="line.127"></a>
131 <FONT color="green">128</FONT> // Cross partial derivatives<a name="line.128"></a>
132 <FONT color="green">129</FONT> final double[][] dZdXdY = new double[xLen][yLen];<a name="line.129"></a>
133 <FONT color="green">130</FONT> for (int i = 0; i &lt; xLen ; i++) {<a name="line.130"></a>
134 <FONT color="green">131</FONT> final int nI = nextIndex(i, xLen);<a name="line.131"></a>
135 <FONT color="green">132</FONT> final int pI = previousIndex(i);<a name="line.132"></a>
136 <FONT color="green">133</FONT> for (int j = 0; j &lt; yLen; j++) {<a name="line.133"></a>
137 <FONT color="green">134</FONT> final int nJ = nextIndex(j, yLen);<a name="line.134"></a>
138 <FONT color="green">135</FONT> final int pJ = previousIndex(j);<a name="line.135"></a>
139 <FONT color="green">136</FONT> dZdXdY[i][j] = (zY_2[nI][nJ] - zY_2[nI][pJ] -<a name="line.136"></a>
140 <FONT color="green">137</FONT> zY_2[pI][nJ] + zY_2[pI][pJ]) /<a name="line.137"></a>
141 <FONT color="green">138</FONT> ((xval[nI] - xval[pI]) * (yval[nJ] - yval[pJ])) ;<a name="line.138"></a>
142 <FONT color="green">139</FONT> }<a name="line.139"></a>
143 <FONT color="green">140</FONT> }<a name="line.140"></a>
144 <FONT color="green">141</FONT> <a name="line.141"></a>
145 <FONT color="green">142</FONT> // Create the interpolating splines<a name="line.142"></a>
146 <FONT color="green">143</FONT> return new BicubicSplineInterpolatingFunction(xval, yval, zY_2,<a name="line.143"></a>
147 <FONT color="green">144</FONT> dZdX, dZdY, dZdXdY);<a name="line.144"></a>
148 <FONT color="green">145</FONT> }<a name="line.145"></a>
149 <FONT color="green">146</FONT> <a name="line.146"></a>
150 <FONT color="green">147</FONT> /**<a name="line.147"></a>
151 <FONT color="green">148</FONT> * Compute the next index of an array, clipping if necessary.<a name="line.148"></a>
152 <FONT color="green">149</FONT> * It is assumed (but not checked) that {@code i} is larger than or equal to 0}.<a name="line.149"></a>
153 <FONT color="green">150</FONT> *<a name="line.150"></a>
154 <FONT color="green">151</FONT> * @param i Index<a name="line.151"></a>
155 <FONT color="green">152</FONT> * @param max Upper limit of the array<a name="line.152"></a>
156 <FONT color="green">153</FONT> * @return the next index<a name="line.153"></a>
157 <FONT color="green">154</FONT> */<a name="line.154"></a>
158 <FONT color="green">155</FONT> private int nextIndex(int i, int max) {<a name="line.155"></a>
159 <FONT color="green">156</FONT> final int index = i + 1;<a name="line.156"></a>
160 <FONT color="green">157</FONT> return index &lt; max ? index : index - 1;<a name="line.157"></a>
161 <FONT color="green">158</FONT> }<a name="line.158"></a>
162 <FONT color="green">159</FONT> /**<a name="line.159"></a>
163 <FONT color="green">160</FONT> * Compute the previous index of an array, clipping if necessary.<a name="line.160"></a>
164 <FONT color="green">161</FONT> * It is assumed (but not checked) that {@code i} is smaller than the size of the array.<a name="line.161"></a>
165 <FONT color="green">162</FONT> *<a name="line.162"></a>
166 <FONT color="green">163</FONT> * @param i Index<a name="line.163"></a>
167 <FONT color="green">164</FONT> * @return the previous index<a name="line.164"></a>
168 <FONT color="green">165</FONT> */<a name="line.165"></a>
169 <FONT color="green">166</FONT> private int previousIndex(int i) {<a name="line.166"></a>
170 <FONT color="green">167</FONT> final int index = i - 1;<a name="line.167"></a>
171 <FONT color="green">168</FONT> return index &gt;= 0 ? index : 0;<a name="line.168"></a>
172 <FONT color="green">169</FONT> }<a name="line.169"></a>
173 <FONT color="green">170</FONT> }<a name="line.170"></a>
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234 </PRE>
235 </BODY>
236 </HTML>