> this bug's URL: http://intevation.de/rt/webrt?serial_num=3166
> -------------------------------------------------------------------
>
> Subject: bad target region calculations with i.rectify order=3
[1st and 2nd order work fine]
>
> i.rectify (both 6.0.0, 6.1cvs) is creating rasters with incorrect
> regions and resolutions when doing a 3rd order transform. It works
> fine with 1st or second order transforms. 'i.rectify -c' (or lower
> order) is needed to get a good transform.
>
> Starting with unreferenced simple XY map (PNG imported with
> r.in.gdal) and re-regioned to positive values with r.region, keeping
> res=1 pixel.
>
> Target location is lat/lon (both wgs84 and none/intl datums fail)
>
> e.g. eastern boundary ends up being +2 degrees east of what it
> should be, resolution is calculated from e-w/cols so that ends up
> very coarse... southern boundary is wrong too.
>
> CRS_georef() in imagery/i.rectify/crs.c ?
>
> more info available upon request.
Hamish,
you may compare with GDAL:
gdal/alg/gdal_crs.c
ok.
- is the algorithm identical?
I think so.
The CRS_georef() fns are essentially the same (GRASS 6/GDAL 1.2.5).
I have cut the top off of both files to CRS_compute_georef_equations(),
diff -wi attached.. No difference there that I can see.
Note that Frank has some comments in his:
* $Log: gdal_crs.c,v $
[gdal 1.2.6]
+ * Revision 1.8 2004/12/26 16:12:21 fwarmerdam
+ * thin plate spline support now implemented
+ *
[nice..]
* Revision 1.7 2004/08/11 19:01:56 warmerda
* default to 2nd order if we have enough points for 3rd order but 0 requested
*
* Revision 1.6 2004/08/09 14:38:27 warmerda
* added serialize/deserialize support for warpoptions and transformers
*
* Revision 1.5 2004/03/19 15:22:02 warmerda
* Fixed double free of padfRasterX. Submitted by Scott Reynolds.
*
* Revision 1.4 2002/12/07 17:09:50 warmerda
* re-enable 3rd order even though unstable
*
* Revision 1.3 2002/12/06 17:58:00 warmerda
* fix a few bugs
*
* Revision 1.2 2002/12/05 05:46:08 warmerda
* fixed linkage problem
*
* Revision 1.1 2002/12/05 05:43:23 warmerda
* New
which tells us 3rd order is not stable...?
If we can get the inputs out, maybe we can try doing the matrix math
in Matlab/Octave and see what answers it gives... (polyfit() fn.)
- does gdalwarp produce the same (error)
I am starting with a flat PNG(Tiff) image with no georefencing.
I assume you can use something in GDAL to insert GCPs as defined
in the GRASS group/POINTS file, but not the four corners of the
map - that's where it all goes wrong in i.rectify.
So I don't know how to try gdalwarp on it.
?
If GDAL performs better for above operation, its improvements
might be worth backporting.
I think GDAL will have the same.
?
Hamish
(attachments)
crs.diff (1.53 KB)