svn-gvsig-desktop / tags / v10_RC2c / extensions / extGeoreferencing / src / org / gvsig / georeferencing / GeoOperations.java @ 8745
History | View | Annotate | Download (20.6 KB)
1 | 5791 | nacho | /* gvSIG. Sistema de Informaci?n Geogr?fica de la Generalitat Valenciana
|
---|---|---|---|
2 | *
|
||
3 | * Copyright (C) 2006 IVER T.I. and Generalitat Valenciana.
|
||
4 | *
|
||
5 | * This program is free software; you can redistribute it and/or
|
||
6 | * modify it under the terms of the GNU General Public License
|
||
7 | * as published by the Free Software Foundation; either version 2
|
||
8 | * of the License, or (at your option) any later version.
|
||
9 | *
|
||
10 | * This program is distributed in the hope that it will be useful,
|
||
11 | * but WITHOUT ANY WARRANTY; without even the implied warranty of
|
||
12 | * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
||
13 | * GNU General Public License for more details.
|
||
14 | *
|
||
15 | * You should have received a copy of the GNU General Public License
|
||
16 | * along with this program; if not, write to the Free Software
|
||
17 | * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,USA.
|
||
18 | */
|
||
19 | 5217 | nacho | package org.gvsig.georeferencing; |
20 | |||
21 | import java.awt.geom.Point2D; |
||
22 | import java.io.BufferedOutputStream; |
||
23 | 5328 | nacho | import java.io.BufferedReader; |
24 | import java.io.BufferedWriter; |
||
25 | 5217 | nacho | import java.io.DataOutputStream; |
26 | import java.io.File; |
||
27 | 5328 | nacho | import java.io.FileInputStream; |
28 | import java.io.FileNotFoundException; |
||
29 | 5217 | nacho | import java.io.FileOutputStream; |
30 | import java.io.IOException; |
||
31 | 5328 | nacho | import java.io.InputStreamReader; |
32 | import java.io.OutputStream; |
||
33 | import java.io.OutputStreamWriter; |
||
34 | 5217 | nacho | import java.util.zip.DataFormatException; |
35 | |||
36 | 5982 | nacho | import org.cresques.io.data.RasterMetaFileTags; |
37 | 5217 | nacho | import org.gvsig.georeferencing.utils.AffineT; |
38 | 5697 | nacho | import org.gvsig.georeferencing.utils.GeoUtils; |
39 | 5217 | nacho | import org.gvsig.georeferencing.utils.PointT; |
40 | import org.xmlpull.v1.XmlPullParserException; |
||
41 | import org.xmlpull.v1.XmlPullParserFactory; |
||
42 | import org.xmlpull.v1.XmlSerializer; |
||
43 | |||
44 | 8546 | nacho | import com.iver.cit.gvsig.fmap.DriverException; |
45 | 5241 | nacho | import com.iver.cit.gvsig.fmap.layers.FLyrPoints; |
46 | 5217 | nacho | |
47 | 5241 | nacho | |
48 | 5217 | nacho | /**
|
49 | * Operaciones necesarias para la georreferenciaci?n a partir de la capa de puntos.
|
||
50 | * En base a unos puntos de control de origen y otros de destino se crea una matriz
|
||
51 | * de transformaci?n para asignar la nueva posici?n a la imagen y crear ficheros
|
||
52 | * de georreferenciaci?n en dos formatos: worldfile y rasterMetaFile
|
||
53 | *
|
||
54 | * @author Nacho Brodin (brodin_ign@gva.es)
|
||
55 | */
|
||
56 | public class GeoOperations{ |
||
57 | |||
58 | //**********************Vars**********************************
|
||
59 | private int order = 1; |
||
60 | private PointT[] src = null; |
||
61 | private PointT[] dst = null; |
||
62 | 8546 | nacho | private AffineT affine = null; |
63 | private boolean createWorldFile = false; |
||
64 | 5217 | nacho | //**********************End Vars******************************
|
65 | |||
66 | //**********************Methods*******************************
|
||
67 | /**
|
||
68 | * Constructor
|
||
69 | */
|
||
70 | public GeoOperations(){}
|
||
71 | |||
72 | /**
|
||
73 | * Constructor. Crea las estructuras de puntos y calcula la transformaci?n af?n.
|
||
74 | * @param lyr
|
||
75 | */
|
||
76 | public GeoOperations(FLyrPoints lyr){
|
||
77 | FLyrPoints lyrPoints = lyr; |
||
78 | 5491 | nacho | src = new PointT[lyr.getCountActivePoints()];
|
79 | dst = new PointT[lyr.getCountActivePoints()];
|
||
80 | 5217 | nacho | int nPoint = 0; |
81 | 6097 | nacho | for(int i = 0; i<lyr.getCountPoints(); i++){ |
82 | 5217 | nacho | if(lyr.getPoint(i).active == true){ |
83 | src[nPoint] = new PointT();
|
||
84 | dst[nPoint] = new PointT();
|
||
85 | src[nPoint].setX(lyr.getPoint(i).pixelPoint.getX()); |
||
86 | src[nPoint].setY(lyr.getPoint(i).pixelPoint.getY()); |
||
87 | src[nPoint].setI(lyr.getCountPoints()); |
||
88 | dst[nPoint].setX(lyr.getPoint(i).mapPoint.getX()); |
||
89 | dst[nPoint].setY(lyr.getPoint(i).mapPoint.getY()); |
||
90 | dst[nPoint].setI(lyr.getCountPoints()); |
||
91 | nPoint++; |
||
92 | } |
||
93 | } |
||
94 | 8546 | nacho | |
95 | 5491 | nacho | if(lyr.getCountActivePoints() >= 3 * 10) |
96 | 5217 | nacho | order = 3;
|
97 | 5491 | nacho | else if(lyr.getCountActivePoints() >= 3 * 6) |
98 | 5217 | nacho | order = 2;
|
99 | else
|
||
100 | order = 1;
|
||
101 | |||
102 | affine = new AffineT();
|
||
103 | |||
104 | //Calcular la transformaci?n afin por m?nimos cuadrados
|
||
105 | affine = computeLeastSquaresAffine(affine, src, dst, order); |
||
106 | } |
||
107 | |||
108 | /**
|
||
109 | * A partir de la transformaci?n af?n creada en el construtor
|
||
110 | * genera un fichero de georreferenciaci?n para la imagen.
|
||
111 | * @param lyr Capa de puntos
|
||
112 | * @param order Orden del sistema de ecuaciones
|
||
113 | * @param widthPx Ancho en pixeles de la imagen a georreferenciar
|
||
114 | * @param heightPx Alto en pixeles de la imagen a georreferenciar
|
||
115 | * @param file Nombre del fichero raster a georreferenciar
|
||
116 | */
|
||
117 | public void createGeorefFile(int widthPx, int heightPx, String file){ |
||
118 | try{
|
||
119 | File f = new File(file); |
||
120 | 6214 | nacho | String nameWorldFile = f.getPath().substring(0, f.getPath().lastIndexOf(".")) + GeoUtils.getWorldFileExtensionFromFileName(file); |
121 | 5217 | nacho | if(createWorldFile)
|
122 | createWorldFile(affine, widthPx, heightPx, nameWorldFile); |
||
123 | createRasterMetaFile(affine, widthPx, heightPx, nameWorldFile.substring(0, nameWorldFile.lastIndexOf("."))+".rmf"); |
||
124 | }catch(IOException ex){ |
||
125 | System.err.println("Can't create WorldFile"); |
||
126 | } |
||
127 | } |
||
128 | |||
129 | /**
|
||
130 | * A partir de la transformaci?n af?n calculada en el contructor
|
||
131 | * transforma los puntos pasados como par?metros.
|
||
132 | * @param lyr Capa de puntos
|
||
133 | * @param list Lista de puntos a transformar
|
||
134 | * @return Lista de puntos transformados
|
||
135 | */
|
||
136 | public Point2D[] transformPoints(Point2D[] list){ |
||
137 | Point2D[] result = new Point2D[list.length]; |
||
138 | for(int i = 0; i < list.length;i++){ |
||
139 | double[] pto = transformPoint((int)list[i].getX(), (int)list[i].getY(), affine); |
||
140 | result[i] = new Point2D.Double(pto[0], pto[1]); |
||
141 | } |
||
142 | return result;
|
||
143 | } |
||
144 | |||
145 | /**
|
||
146 | * Crea un fichero de georreferenciaci?n (worldfile) a partir de la transformaci?n af?n. Para
|
||
147 | * esto necesita obtener las coordenadas reales de la coordenada en pixels (0,0) y calcular
|
||
148 | * el tama?o de pixel.
|
||
149 | * @param affine Transformaci?n
|
||
150 | * @param widthPx Ancho en pixeles de la imagen a georreferenciar
|
||
151 | * @param heightPx Alto en pixeles de la imagen a georreferenciar
|
||
152 | * @param nameWordlFile Nombre del fichero de georreferenciaci?n
|
||
153 | * @throws IOException
|
||
154 | */
|
||
155 | private void createWorldFile(AffineT affine, int widthPx, int heightPx, String nameWordlFile) throws IOException { |
||
156 | StringBuffer data = new StringBuffer(); |
||
157 | 5896 | luisw2 | //double[] begin = transformPoint(0, 0, affine);
|
158 | //double[] end = transformPoint(widthPx, heightPx, affine);
|
||
159 | //double pixelSizeX = (end[0] - begin[0])/(widthPx - 1);
|
||
160 | //double pixelSizeY = (end[1] - begin[1])/(heightPx - 1);
|
||
161 | |||
162 | 5217 | nacho | File f = new File(nameWordlFile); |
163 | DataOutputStream dos = new DataOutputStream( new BufferedOutputStream(new FileOutputStream(f)) ); |
||
164 | |||
165 | 5896 | luisw2 | data.append(affine.getCofX(1)+"\n");//pixelSizeX+"\n"); |
166 | data.append(affine.getCofX(2)+"\n"); |
||
167 | data.append(affine.getCofY(1)+"\n"); |
||
168 | data.append(affine.getCofY(2)+"\n");//(-1 * pixelSizeY)+"\n"); |
||
169 | data.append(affine.getCofX(0)+"\n");//""+begin[0]+"\n"); |
||
170 | data.append(affine.getCofY(0)+"\n");//""+end[1]+"\n"); |
||
171 | 5217 | nacho | |
172 | dos.writeBytes(data.toString()); |
||
173 | dos.close(); |
||
174 | } |
||
175 | |||
176 | /**
|
||
177 | 5328 | nacho | * A partir de XmlSerializer se salva la georreferenciaci?n
|
178 | * @param serializer
|
||
179 | * @param min
|
||
180 | * @param max
|
||
181 | * @param widthPx Ancho en pixeles de la imagen a georreferenciar
|
||
182 | * @param heightPx Alto en pixeles de la imagen a georreferenciar
|
||
183 | * @throws IOException
|
||
184 | */
|
||
185 | private void serializeGeoreferencing(XmlSerializer serializer, double[] min, double[] max, int widthPx, int heightPx) throws IOException { |
||
186 | double pixelSizeX = (max[0] - min[0])/(widthPx - 1); |
||
187 | double pixelSizeY = (max[1] - min[1])/(heightPx - 1); |
||
188 | |||
189 | serializer.startTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.PROJ).text("Projection").endTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.PROJ).text("\n"); |
||
190 | serializer.startTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.BBOX).text("\n");
|
||
191 | serializer.startTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.POSX).text(""+min[0]).endTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.POSX).text("\n"); |
||
192 | 5952 | nacho | serializer.startTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.POSY).text(""+min[1]).endTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.POSY).text("\n"); |
193 | 5982 | nacho | serializer.startTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.ROTX).text(""+affine.getCofX(2)).endTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.ROTX).text("\n"); |
194 | serializer.startTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.ROTY).text(""+affine.getCofY(1)).endTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.ROTY).text("\n"); |
||
195 | 5328 | nacho | serializer.startTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.PX_SIZE_X).text(""+pixelSizeX).endTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.PX_SIZE_X).text("\n"); |
196 | 5952 | nacho | serializer.startTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.PX_SIZE_Y).text(""+ pixelSizeY).endTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.PX_SIZE_Y).text("\n"); |
197 | 5328 | nacho | serializer.startTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.WIDTH).text(""+(max[0] - min[0])).endTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.WIDTH).text("\n"); |
198 | serializer.startTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.HEIGHT).text(""+(max[1] - min[1])).endTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.HEIGHT).text("\n"); |
||
199 | serializer.endTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.BBOX).text("\n");
|
||
200 | serializer.startTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.DIM).text("\n");
|
||
201 | serializer.startTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.PX_WIDTH).text(""+widthPx).endTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.PX_WIDTH).text("\n"); |
||
202 | serializer.startTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.PX_HEIGHT).text(""+heightPx).endTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.PX_HEIGHT).text("\n"); |
||
203 | serializer.endTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.DIM).text("\n");
|
||
204 | } |
||
205 | |||
206 | /**
|
||
207 | * Si ya existe un fichero de georreferenciaci?n se ha creado un fichero .tmp con
|
||
208 | * la nueva georreferenciaci?n. Esta funci?n mezcla este temporal con el .rmf existente
|
||
209 | * en un nuevo fichero _tmpOutput que ser? renombrado como el nuevo .rmf. Finalmente
|
||
210 | * elimina el temporal y el _tmpOutput.
|
||
211 | * @param fileName
|
||
212 | */
|
||
213 | private void mergeFiles(String fileName){ |
||
214 | try{
|
||
215 | File rmfFile = new File(fileName); |
||
216 | File tmpFile = new File(fileName.substring(0, fileName.lastIndexOf("."))+".tmp"); |
||
217 | File out = new File(rmfFile+"_tmpOutput"); |
||
218 | |||
219 | BufferedReader inRmf = new BufferedReader(new InputStreamReader(new FileInputStream(rmfFile))); |
||
220 | BufferedReader inTmp = new BufferedReader(new InputStreamReader(new FileInputStream(tmpFile))); |
||
221 | BufferedWriter outTmp = new BufferedWriter(new OutputStreamWriter(new FileOutputStream(out))); |
||
222 | |||
223 | //Leemos el principio del .rmf
|
||
224 | String str = inRmf.readLine();
|
||
225 | 5352 | nacho | while(!str.startsWith("<"+RasterMetaFileTags.MAIN_TAG)){ |
226 | 5328 | nacho | outTmp.write(str+"\n");
|
227 | str = inRmf.readLine(); |
||
228 | } |
||
229 | 5352 | nacho | outTmp.write(str+"\n");
|
230 | 5328 | nacho | |
231 | //Leemos la georreferenciaci?n del .tmp
|
||
232 | 5352 | nacho | while(str != null && !str.startsWith("<"+RasterMetaFileTags.LAYER)) |
233 | 5328 | nacho | str = inTmp.readLine(); |
234 | while(str != null){ |
||
235 | outTmp.write(str+"\n");
|
||
236 | str = inTmp.readLine(); |
||
237 | } |
||
238 | |||
239 | //Saltamos la georreferenciaci?n que ya existia en el .rmf
|
||
240 | 5352 | nacho | try{
|
241 | 5328 | nacho | str = inRmf.readLine(); |
242 | 5352 | nacho | while(str != null && |
243 | !str.startsWith("</"+RasterMetaFileTags.LAYER) &&
|
||
244 | !str.startsWith("<"+RasterMetaFileTags.GEOPOINTS) &&
|
||
245 | !str.startsWith("</"+RasterMetaFileTags.MAIN_TAG))
|
||
246 | str = inRmf.readLine(); |
||
247 | }catch(Exception exc){ |
||
248 | |||
249 | } |
||
250 | 5328 | nacho | |
251 | //Leemos el resto del fichero .rmf
|
||
252 | 5352 | nacho | if(!str.startsWith("<"+RasterMetaFileTags.GEOPOINTS)) |
253 | str = inRmf.readLine(); |
||
254 | 5328 | nacho | while(str != null){ |
255 | outTmp.write(str+"\n");
|
||
256 | str = inRmf.readLine(); |
||
257 | } |
||
258 | |||
259 | outTmp.close(); |
||
260 | inRmf.close(); |
||
261 | inTmp.close(); |
||
262 | |||
263 | //Eliminamos el antiguo .rmf y lo sustituimos
|
||
264 | rmfFile.delete(); |
||
265 | out.renameTo(rmfFile); |
||
266 | tmpFile.delete(); |
||
267 | |||
268 | }catch(FileNotFoundException exc){ |
||
269 | |||
270 | }catch(IOException exc){ |
||
271 | |||
272 | } |
||
273 | } |
||
274 | |||
275 | /**
|
||
276 | 5217 | nacho | * Si no existe crea un fichero .rmf con la georreferenciaci?n de la imagen. Si existe este
|
277 | * deber? a?adir o modificar la informaci?n de georreferenciaci?n en el fichero.
|
||
278 | * Para esto necesita obtener las coordenadas reales de la coordenada en pixels (0,0) y
|
||
279 | * calcular el tama?o de pixel.
|
||
280 | * @param affine Transformaci?n
|
||
281 | * @param widthPx Ancho en pixeles de la imagen a georreferenciar
|
||
282 | * @param heightPx Alto en pixeles de la imagen a georreferenciar
|
||
283 | * @param nameWordlFile Nombre del fichero de georreferenciaci?n
|
||
284 | * @throws IOException
|
||
285 | */
|
||
286 | 5328 | nacho | private void createRasterMetaFile(AffineT affine, int widthPx, int heightPx, String nameWorldFile) throws IOException { |
287 | File file = new File(nameWorldFile); |
||
288 | 5217 | nacho | double[] min = transformPoint(0, 0, affine); |
289 | double[] max = transformPoint(widthPx, heightPx, affine); |
||
290 | 5328 | nacho | try{
|
291 | XmlPullParserFactory factory = XmlPullParserFactory.newInstance(System.getProperty(XmlPullParserFactory.PROPERTY_NAME), null); |
||
292 | XmlSerializer serializer = factory.newSerializer(); |
||
293 | 5217 | nacho | |
294 | 5328 | nacho | if(file.exists()){
|
295 | file = new File(nameWorldFile.substring(0, nameWorldFile.lastIndexOf("."))+".tmp"); |
||
296 | 5217 | nacho | serializer.setOutput(new FileOutputStream(file), null); |
297 | serializer.startDocument(null, null); |
||
298 | serializer.ignorableWhitespace("\n\n");
|
||
299 | serializer.setPrefix("", RasterMetaFileTags.NAMESPACE);
|
||
300 | serializer.startTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.LAYER).text("\n");
|
||
301 | 5328 | nacho | serializeGeoreferencing(serializer, min, max, widthPx, heightPx); |
302 | 5217 | nacho | serializer.endTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.LAYER).text("\n");
|
303 | serializer.endDocument(); |
||
304 | 5328 | nacho | mergeFiles(nameWorldFile); |
305 | }else{
|
||
306 | serializer.setOutput(new FileOutputStream(file), null); |
||
307 | serializer.startDocument(null, null); |
||
308 | serializer.ignorableWhitespace("\n\n");
|
||
309 | serializer.setPrefix("", RasterMetaFileTags.NAMESPACE);
|
||
310 | serializer.startTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.MAIN_TAG).text("\n");
|
||
311 | serializer.startTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.LAYER).text("\n");
|
||
312 | serializeGeoreferencing(serializer, min, max, widthPx, heightPx); |
||
313 | serializer.endTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.LAYER).text("\n");
|
||
314 | serializer.endTag(RasterMetaFileTags.NAMESPACE, RasterMetaFileTags.MAIN_TAG).text("\n");
|
||
315 | serializer.endDocument(); |
||
316 | } |
||
317 | }catch(XmlPullParserException exc){}
|
||
318 | 5217 | nacho | } |
319 | |||
320 | /**
|
||
321 | * Calcula la transformaci?n af?n a partir de los puntos introducidos por el m?todo de
|
||
322 | * m?nimos cuadrados.
|
||
323 | * @param af Transformaci?n af?n
|
||
324 | * @param src Puntos de entrada en coordenadas pixel
|
||
325 | * @param dst Puntos de destino en coordenadas del mundo
|
||
326 | * @param order Orden del sistema de ecuaciones
|
||
327 | */
|
||
328 | public AffineT computeLeastSquaresAffine(AffineT af, PointT[] src, PointT[] dst, int order){ |
||
329 | double[] b = new double[dst.length]; |
||
330 | int i,cofs;
|
||
331 | |||
332 | af.setOrder(order); |
||
333 | cofs = order <= 2 ? order*3 : 10; |
||
334 | af.mallocCofs(cofs); |
||
335 | |||
336 | 8546 | nacho | if(src.length == 2){ |
337 | af.setCofX(0, 2); |
||
338 | af.setCofY(0, 1); |
||
339 | double distWCX = dst[1].getX() - dst[0].getX(); |
||
340 | double distWCY = dst[1].getY() - dst[0].getY(); |
||
341 | double distPxX = src[1].getX() - src[0].getX(); |
||
342 | double distPxY = src[1].getY() - src[0].getY(); |
||
343 | double psX = distWCX / distPxX;
|
||
344 | double psY = distWCY / distPxY;
|
||
345 | af.setCofX(psX, 1);
|
||
346 | af.setCofY(psY, 2);
|
||
347 | |||
348 | double initWCX = dst[0].getX() - ((src[0].getX() + 1) * psX); |
||
349 | double initWCY = dst[0].getY() - ((src[0].getY() + 1) * psY); |
||
350 | |||
351 | af.setCofX(initWCX, 0);
|
||
352 | af.setCofY(initWCY, 0);
|
||
353 | return af;
|
||
354 | } |
||
355 | |||
356 | 5217 | nacho | //First compute the X cofs
|
357 | for(i = 0; i < dst.length; i++) { |
||
358 | b[i] = dst[i].getX(); |
||
359 | } |
||
360 | |||
361 | af.setCofX(singleLeastSquaresAffine(af.getXcofs(), src, b, order)); |
||
362 | |||
363 | //Now compute the Y cofs
|
||
364 | for(i = 0; i < dst.length; i++) { |
||
365 | b[i] = dst[i].getY(); |
||
366 | } |
||
367 | |||
368 | af.setCofY(singleLeastSquaresAffine(af.getYcofs(), src, b, order)); |
||
369 | |||
370 | return af;
|
||
371 | } |
||
372 | |||
373 | private double[] singleLeastSquaresAffine(double[] a, PointT[] src, double[] dst, int order){ |
||
374 | int i,j;
|
||
375 | int n = dst.length;
|
||
376 | |||
377 | int points = order <= 2 ? order * 3 : 10; |
||
378 | double[][] combined = new double[points][points + 1]; |
||
379 | double[][] A = new double[n + 1][points]; |
||
380 | double[][] b = new double[n + 1][1]; |
||
381 | |||
382 | for(i = 0; i < n; i++) { |
||
383 | A[i][0] = 1.0; |
||
384 | A[i][1] = src[i].getX();
|
||
385 | A[i][2] = src[i].getY();
|
||
386 | if(order > 1) { |
||
387 | A[i][3] = src[i].getX() * src[i].getX();
|
||
388 | A[i][4] = src[i].getX() * src[i].getY();
|
||
389 | A[i][5] = src[i].getY() * src[i].getY();
|
||
390 | } |
||
391 | if(order > 2) { |
||
392 | A[i][6] = src[i].getX() * src[i].getX() * src[i].getX();
|
||
393 | A[i][7] = src[i].getX() * src[i].getX() * src[i].getY();
|
||
394 | A[i][8] = src[i].getX() * src[i].getY() * src[i].getY();
|
||
395 | A[i][9] = src[i].getY() * src[i].getY() * src[i].getY();
|
||
396 | } |
||
397 | b[i][0] = dst[i];
|
||
398 | } |
||
399 | |||
400 | 5697 | nacho | double[][] Atrans = GeoUtils.transpose(A, n, points); |
401 | double[][] left = GeoUtils.multmatrix(Atrans,A,points,n,points); |
||
402 | double[][] right = GeoUtils.multmatrix(Atrans,b,points,n,1); |
||
403 | 5217 | nacho | |
404 | for(i = 0; i < points; i++) { |
||
405 | combined[i][0] = right[i][0]; |
||
406 | for(j = 0; j < points; j++) { |
||
407 | combined[i][j+1] = left[i][j];
|
||
408 | } |
||
409 | } |
||
410 | |||
411 | try{
|
||
412 | combined = solve(combined, points, points+1);
|
||
413 | }catch(DataFormatException ex){ |
||
414 | System.err.println("Can't solve matrix"); |
||
415 | } |
||
416 | |||
417 | for(i = 0; i < points; i++) { |
||
418 | a[i] = combined[i][0];
|
||
419 | } |
||
420 | return a;
|
||
421 | } |
||
422 | |||
423 | |||
424 | private double[][] solve(double[][] mat,int rows,int cols)throws DataFormatException{ |
||
425 | int i,j,k;
|
||
426 | double d,tmp;
|
||
427 | int big;
|
||
428 | |||
429 | for(i = 0; i < rows; i++) { |
||
430 | // Find largest row
|
||
431 | big = i; |
||
432 | for(j = i; j < rows; j++) {
|
||
433 | if(Math.abs(mat[j][i+1]) > Math.abs(mat[big][i+1])) { |
||
434 | big = j; |
||
435 | } |
||
436 | } |
||
437 | // swap row i and row big
|
||
438 | for(k = 0; k < cols ; k++) { |
||
439 | tmp = mat[i][k]; |
||
440 | mat[i][k] = mat[big][k]; |
||
441 | mat[big][k] = tmp; |
||
442 | } |
||
443 | if(mat[i][i+1] == 0) |
||
444 | throw new DataFormatException(); |
||
445 | |||
446 | d = 1.0 / mat[i][i+1]; |
||
447 | for(j = 0; j < cols ; j++) { |
||
448 | mat[i][j] *= d; |
||
449 | //assert(!isnan(mat[i][j]));
|
||
450 | } |
||
451 | for(k = 0; k < rows; k++) { |
||
452 | if(k == i)
|
||
453 | continue;
|
||
454 | if(mat[k][i+1] != 0.0) { |
||
455 | d = mat[k][i+1] / mat[i][i+1]; |
||
456 | for(j = 0; j < cols ; j++) { |
||
457 | mat[k][j] -= d*mat[i][j]; |
||
458 | //assert(!isnan(mat[k][j]));
|
||
459 | } |
||
460 | } |
||
461 | } |
||
462 | } |
||
463 | return mat;
|
||
464 | } |
||
465 | |||
466 | /**
|
||
467 | * Obtiene las coordenadas de un punto en coordenadas del mundo a partir de una transformaci?n
|
||
468 | * y el punto de entrada en coordenadas de la imagen.
|
||
469 | * @param sx Coordenada X del punto de entrada
|
||
470 | * @param syz Coordenada Y del punto de entrada
|
||
471 | * @param af Transformaci?n
|
||
472 | * @return Punto transformado
|
||
473 | */
|
||
474 | public double[] transformPoint(int sx, int sy, AffineT af){ |
||
475 | double[] p = new double[2]; |
||
476 | if(af.getOrder() == 1) { |
||
477 | p[0] = af.getCofX(0) + af.getCofX(1) * sx + af.getCofX(2) * sy; |
||
478 | p[1] = af.getCofY(0) + af.getCofY(1) * sx + af.getCofY(2) * sy; |
||
479 | } |
||
480 | else if(af.getOrder() == 2) { |
||
481 | p = quadTransformPoint(sx, sy, af); |
||
482 | } |
||
483 | else {
|
||
484 | p = cubicTransformPoint(sx, sy, af); |
||
485 | } |
||
486 | return p;
|
||
487 | } |
||
488 | |||
489 | private static double[] quadTransformPoint(int sx, int sy, AffineT af){ |
||
490 | double[] p = new double[2]; |
||
491 | p[0] = af.getCofX(0) + af.getCofX(1) * sx + af.getCofX(2) * sy + |
||
492 | af.getCofX(3) * sx * sx + af.getCofX(4) * sx * sy + af.getCofX(5) * sy * sy; |
||
493 | p[1] = af.getCofY(0) + af.getCofY(1) * sx + af.getCofY(2) * sy + |
||
494 | af.getCofY(3) * sx * sx + af.getCofY(4) * sx * sy + af.getCofY(5) * sy * sy; |
||
495 | return p;
|
||
496 | } |
||
497 | |||
498 | private static double[] cubicTransformPoint(int sx, int sy, AffineT af){ |
||
499 | double[] p = new double[2]; |
||
500 | p[0] = af.getCofX(0) + af.getCofX(1) * sx + af.getCofX(2) * sy + |
||
501 | af.getCofX(3) * sx * sx + af.getCofX(4) * sx * sy + af.getCofX(5) * sy * sy + |
||
502 | af.getCofX(6) * sx * sx * sx + af.getCofX(7) * sx * sx * sy + |
||
503 | af.getCofX(8)*sx*sy*sy + af.getCofX(9)*sy*sy*sy; |
||
504 | p[1] = af.getCofY(0) + af.getCofY(1) * sx + af.getCofY(2) * sy + |
||
505 | af.getCofY(3) * sx * sx + af.getCofY(4) * sx * sy + af.getCofY(5) * sy * sy + |
||
506 | af.getCofY(6) * sx * sx * sx + af.getCofY(7) * sx * sx * sy + |
||
507 | af.getCofY(8) * sx * sy * sy + af.getCofY(9) * sy * sy * sy; |
||
508 | return p;
|
||
509 | } |
||
510 | //**********************End Methods***************************
|
||
511 | |||
512 | //**********************Setters & Getters*********************
|
||
513 | /**
|
||
514 | * M?todo para notificar si se desea crear o no el fichero de coordenadas .tfw, .wld .jpgw ,...
|
||
515 | * @param createWorldFile
|
||
516 | */
|
||
517 | public void setCreateWorldFile(boolean createWorldFile) { |
||
518 | this.createWorldFile = createWorldFile;
|
||
519 | } |
||
520 | |||
521 | /**
|
||
522 | * Obtiene la extensi?n del fichero de georreferenciaci?n
|
||
523 | * @return String con la extensi?n del fichero de georreferenciaci?n dependiendo
|
||
524 | * del valor del formato obtenido del servidor. Por defecto asignaremos un .wld
|
||
525 | */
|
||
526 | 6214 | nacho | /*private String getExtensionWorldFile(String file){
|
527 | 5896 | luisw2 | String ext = file.substring(file.lastIndexOf(".") + 1).toLowerCase();
|
528 | 5217 | nacho | String extWorldFile = ".wld";
|
529 | if(ext.equals("tif") || ext.equals("tiff"))
|
||
530 | extWorldFile = ".tfw";
|
||
531 | if(ext.equals("jpeg") || ext.equals("jpg"))
|
||
532 | 5896 | luisw2 | extWorldFile = ".jgw";
|
533 | if(ext.equals("gif"))
|
||
534 | extWorldFile = ".gfw";
|
||
535 | if(ext.equals("png"))
|
||
536 | extWorldFile = ".pgw";
|
||
537 | 5217 | nacho | return extWorldFile;
|
538 | 6214 | nacho | }*/
|
539 | 5445 | nacho | |
540 | /**
|
||
541 | * Obtiene la matriz de transformaci?n
|
||
542 | * @return AffineT
|
||
543 | */
|
||
544 | public AffineT getAffine() {
|
||
545 | return affine;
|
||
546 | } |
||
547 | 5217 | nacho | //**********************End Setters & Getters*****************
|
548 | 5445 | nacho | |
549 | 8546 | nacho | } |