Quantcast

[gdal-dev] adding georeferencing information

classic Classic list List threaded Threaded
14 messages Options
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

[gdal-dev] adding georeferencing information

Jan Hartmann
Hi all,

Is there a quick way to add back georeferencing information to a geotif file that has been processed by Imagemagick? I have lots of georeferenced images that I want to sharpen up a bit, but ImageMagick throws the georeference headers away. Is there a way to transport them from the original file to the derived one?

Jan

_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: adding georeferencing information

Brent Fraser
Jan,

  I use:
listgeo mygeotiff.tif > mygeotiff.gtt
geotifcp  -g mygeotiff.gtt nogeotags.tif mynewgeotiff.tif


FYI, I've been experimenting with GDAL's vrt format to sharpen my tiffs:

<VRTDataset rasterXSize="2890" rasterYSize="3451">
    :
    :
  <VRTRasterBand dataType="Byte" band="1">
    <Metadata />
    <ColorInterp>Red</ColorInterp>
    <KernelFilteredSource>
      <SourceFilename relativeToVRT="1">NRT_123.tif</SourceFilename>
      <SourceBand>1</SourceBand>
      <SourceProperties RasterXSize="2890" RasterYSize="3451" DataType="Byte" BlockXSize="2890" BlockYSize="1" />
      <SrcRect xOff="0" yOff="0" xSize="2890" ySize="3451" />
      <DstRect xOff="0" yOff="0" xSize="2890" ySize="3451" />
      <Kernel>
        <Size>3</Size>
        <Coefs>-0.111 -0.111 -0.111 -0.111 2 -0.111 -0.111 -0.111 -0.111</Coefs>
      </Kernel>
    </KernelFilteredSource>
  </VRTRasterBand>
  :
  :
</VRTDataset>

My kernel tends to brighten the image a little (note the "2" in the kernel) as well as sharpen.  To try it on your file, create a VRT file with gdal_translate -of VRT, then use a text editor to change the SimpleSource tags to KernelFilteredSource and add the Kernel.  You can then create a sharpened tif by using the VRT as input to gdal_translate, os simply open the VRT in Quantum GIS to view the results.
Best Regards,
Brent Fraser

On 9/5/2011 8:46 AM, Jan Hartmann wrote:
Hi all,

Is there a quick way to add back georeferencing information to a geotif file that has been processed by Imagemagick? I have lots of georeferenced images that I want to sharpen up a bit, but ImageMagick throws the georeference headers away. Is there a way to transport them from the original file to the derived one?

Jan


_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev

_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: adding georeferencing information

Andrew Brooks
In reply to this post by Jan Hartmann
On 5 September 2011 15:46, Jan Hartmann <[hidden email]> wrote:
> Is there a quick way to add back georeferencing information to a geotif file
> that has been processed by Imagemagick?

There should be an "applygeo" utility distributed with the geotiff
library which can do this.

Andrew
_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: adding georeferencing information

Matt Wilkie-2
> Is there a quick way to add back georeferencing information to a geotif file that has been processed by Imagemagick?

also see
http://code.google.com/p/maphew/source/browse/gis/gdal_extras/bin/gdalcopyproj.py

cheers,

matt wilkie
--------------------------------------------
Geomatics Analyst
Information Management and Technology
Yukon Department of Environment
10 Burns Road * Whitehorse, Yukon * Y1A 4Y9
867-667-8133 Tel * 867-393-7003 Fax
http://environmentyukon.gov.yk.ca/geomatics/
--------------------------------------------

_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

[gdal-dev] Numerical instability with thin plate spline transform

Jan Hartmann
In reply to this post by Jan Hartmann
Not sure whether this can be considered a bug, so I give it for what it is worth.

I'm doing thin plate spline transformation from one set of projected coordinates to another. Both sets have values between -600000 and 600000 (meters). A typical set of gcps looks like:

 -gcp 62402 -74383 181915 315371 \
 -gcp -100169 -24213 19685 366661 \
 -gcp 118323 233822 239994 623190 \
...

When transforming the the first two columns of this list with gdaltransform -tps, using the same gcp-list, the results should be exactly equal to the last two columns:

 (from gdal_tps.cpp:)
 * The thin plate spline transformer produces exact transformation
 * at all control points and smoothly varying transformations between
 * control points with greatest influence from local control points.
 * It is suitable for for many applications not well modelled by polynomial
 * transformations.
 *


However, they aren't. With a few control points the differences are only in millimeters, but with a few hundred gcps the differences become meters, and above thousand points the error can get to fifty meters. The problem gets worse when some control points are very near to other ones. I was happy to discover that dividing all numbers by 1000 solves the problem, but there certainly is a numerical instability in the algorithm.

Of course, thin plate spline transformations normally use scan coordinates as input, and these are almost always small numbers. Even so, I can imagine problems with very large rasters and the rubber sheeting applications I am working with.

Jan 



_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: Numerical instability with thin plate spline transform

Big Bear
The linear solver in the TPS routine is naive for any number of
reasons.1,2,3,...  At the same time, one is going to suffer
considerably when the number of control points is in the thousands.
Slow to evaluate so many coefficients for each point.

It is not such a task to improve the routine.  Best to work with
someone who is cognizant per the implementation.  Likely to be
straightforward.  ( Been there, done that. )

Important to know how the implementation is integrated into gdalwarp
before getting too enthusiastic :-/

Thoughts?

1 Looks to be Gauss-Jordan w/o pivoting
2 Never invert a matrix unless required
3 No scaling of the coordinates

On Thu, Sep 8, 2011 at 5:23 AM, Jan Hartmann <[hidden email]> wrote:

> Not sure whether this can be considered a bug, so I give it for what it is
> worth.
>
> I'm doing thin plate spline transformation from one set of projected
> coordinates to another. Both sets have values between -600000 and 600000
> (meters). A typical set of gcps looks like:
>
>  -gcp 62402 -74383 181915 315371 \
>  -gcp -100169 -24213 19685 366661 \
>  -gcp 118323 233822 239994 623190 \
> ...
>
> When transforming the the first two columns of this list with gdaltransform
> -tps, using the same gcp-list, the results should be exactly equal to the
> last two columns:
>
>  (from gdal_tps.cpp:)
>  * The thin plate spline transformer produces exact transformation
>  * at all control points and smoothly varying transformations between
>  * control points with greatest influence from local control points.
>  * It is suitable for for many applications not well modelled by polynomial
>  * transformations.
>  *
>
>
> However, they aren't. With a few control points the differences are only in
> millimeters, but with a few hundred gcps the differences become meters, and
> above thousand points the error can get to fifty meters. The problem gets
> worse when some control points are very near to other ones. I was happy to
> discover that dividing all numbers by 1000 solves the problem, but there
> certainly is a numerical instability in the algorithm.
>
> Of course, thin plate spline transformations normally use scan coordinates
> as input, and these are almost always small numbers. Even so, I can imagine
> problems with very large rasters and the rubber sheeting applications I am
> working with.
>
> Jan
>
>
>
> _______________________________________________
> gdal-dev mailing list
> [hidden email]
> http://lists.osgeo.org/mailman/listinfo/gdal-dev
>
_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: Numerical instability with thin plate spline transform

Even Rouault
Le vendredi 09 septembre 2011 03:23:29, Big Bear a écrit :

> The linear solver in the TPS routine is naive for any number of
> reasons.1,2,3,...  At the same time, one is going to suffer
> considerably when the number of control points is in the thousands.
> Slow to evaluate so many coefficients for each point.
>
> It is not such a task to improve the routine.  Best to work with
> someone who is cognizant per the implementation.  Likely to be
> straightforward.  ( Been there, done that. )
>
> Important to know how the implementation is integrated into gdalwarp
> before getting too enthusiastic :-/
>
> Thoughts?
>
> 1 Looks to be Gauss-Jordan w/o pivoting
> 2 Never invert a matrix unless required
> 3 No scaling of the coordinates

You might try latest GDAL trunk, in particular
http://trac.osgeo.org/gdal/changeset/22876 .

It offers the option to use libarmadillo to speed-up matrix inversion in
thinplatespline.cpp (speed-up when TPS are in the thousands), but perhaps as a
side effect, you'd get also more numerical stability.

>
> On Thu, Sep 8, 2011 at 5:23 AM, Jan Hartmann <[hidden email]> wrote:
> > Not sure whether this can be considered a bug, so I give it for what it
> > is worth.
> >
> > I'm doing thin plate spline transformation from one set of projected
> > coordinates to another. Both sets have values between -600000 and 600000
> > (meters). A typical set of gcps looks like:
> >
> >  -gcp 62402 -74383 181915 315371 \
> >  -gcp -100169 -24213 19685 366661 \
> >  -gcp 118323 233822 239994 623190 \
> > ...
> >
> > When transforming the the first two columns of this list with
> > gdaltransform -tps, using the same gcp-list, the results should be
> > exactly equal to the last two columns:
> >
> >  (from gdal_tps.cpp:)
> >  * The thin plate spline transformer produces exact transformation
> >  * at all control points and smoothly varying transformations between
> >  * control points with greatest influence from local control points.
> >  * It is suitable for for many applications not well modelled by
> > polynomial * transformations.
> >  *
> >
> >
> > However, they aren't. With a few control points the differences are only
> > in millimeters, but with a few hundred gcps the differences become
> > meters, and above thousand points the error can get to fifty meters. The
> > problem gets worse when some control points are very near to other ones.
> > I was happy to discover that dividing all numbers by 1000 solves the
> > problem, but there certainly is a numerical instability in the
> > algorithm.
> >
> > Of course, thin plate spline transformations normally use scan
> > coordinates as input, and these are almost always small numbers. Even
> > so, I can imagine problems with very large rasters and the rubber
> > sheeting applications I am working with.
> >
> > Jan
> >
> >
> >
> > _______________________________________________
> > gdal-dev mailing list
> > [hidden email]
> > http://lists.osgeo.org/mailman/listinfo/gdal-dev
>
> _______________________________________________
> gdal-dev mailing list
> [hidden email]
> http://lists.osgeo.org/mailman/listinfo/gdal-dev
_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: Numerical instability with thin plate spline transform

Jan Hartmann


On 09/09/2011 08:11 AM, Even Rouault wrote:
You might try latest GDAL trunk, in particular http://trac.osgeo.org/gdal/changeset/22876 . It offers the option to use libarmadillo to speed-up matrix inversion in thinplatespline.cpp (speed-up when TPS are in the thousands), but perhaps as a side effect, you'd get also more numerical stability.

I tried to compile the SVN version a few times the last couple of days, but it always broke down: Can someone do a "make clean" in their sandbox and repair the error?

 g++ -fPIC -L/home/a907hart/local/lib -L/usr/lib  gdaldem.o commonutils.o  -L/home/a907hart/src/gdal -lgdal  -larmadillo -L/home/a907hart/local/lib -lgeos_c -L/home/a907hart/local/pgsql/lib -lpq -lpthread -lm -lrt -ldl  -L/home/a907hart/local/lib   -lcurl            -o gdaldem
gdaldem.o:(.data.rel.ro._ZTV21GDALGeneric3x3Dataset[vtable for GDALGeneric3x3Dataset]+0x30): undefined reference to `GDALDataset::CloseDependentDatasets()'
gdaldem.o:(.data.rel.ro._ZTV22GDALColorReliefDataset[vtable for GDALColorReliefDataset]+0x30): undefined reference to `GDALDataset::CloseDependentDatasets()'


Jan

On Thu, Sep 8, 2011 at 5:23 AM, Jan Hartmann [hidden email] wrote:
Not sure whether this can be considered a bug, so I give it for what it
is worth.

I'm doing thin plate spline transformation from one set of projected
coordinates to another. Both sets have values between -600000 and 600000
(meters). A typical set of gcps looks like:

 -gcp 62402 -74383 181915 315371 \
 -gcp -100169 -24213 19685 366661 \
 -gcp 118323 233822 239994 623190 \
...

When transforming the the first two columns of this list with
gdaltransform -tps, using the same gcp-list, the results should be
exactly equal to the last two columns:

 (from gdal_tps.cpp:)
 * The thin plate spline transformer produces exact transformation
 * at all control points and smoothly varying transformations between
 * control points with greatest influence from local control points.
 * It is suitable for for many applications not well modelled by
polynomial * transformations.
 *


However, they aren't. With a few control points the differences are only
in millimeters, but with a few hundred gcps the differences become
meters, and above thousand points the error can get to fifty meters. The
problem gets worse when some control points are very near to other ones.
I was happy to discover that dividing all numbers by 1000 solves the
problem, but there certainly is a numerical instability in the
algorithm.

Of course, thin plate spline transformations normally use scan
coordinates as input, and these are almost always small numbers. Even
so, I can imagine problems with very large rasters and the rubber
sheeting applications I am working with.

Jan



_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev
_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev

    

_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: Numerical instability with thin plate spline transform

multiquadric
As I have indicated previously, a change in method beyond a more
robust linear solver is required.  I have taken the time to fool with
Armadillo.  At the very least I can pass along a trick to improve the
numerical accuracy.. some.

It would help me if you would share the data because I could actually
generate code and an example with real data.  At the same time, I have
never used the TPS option with gdalwarp so it would speed the learning
curve for me (though it looks simple enough).

If you are uncomfortable sharing all the data, the considering sending
the points.  Even as an attachment.  I can do a check to see if you
have a chance with Armadillo w/o additional work.

Please understand that I am not trying to be difficult, but rather
efficient.  I will share comment with Even R. and Frank W. later
Monday or Tuesday with enough information for some triage.  A more
robust solution will follow.  I have the code for this problem that
merely needs some additional notes so that we can all move ahead in a
productive manner.

Best
= David

On Sat, Sep 17, 2011 at 4:10 AM, Jan Hartmann <[hidden email]> wrote:

>
>
> On 09/09/2011 08:11 AM, Even Rouault wrote:
>
> You might try latest GDAL trunk, in particular
> http://trac.osgeo.org/gdal/changeset/22876 . It offers the option to use
> libarmadillo to speed-up matrix inversion in thinplatespline.cpp (speed-up
> when TPS are in the thousands), but perhaps as a side effect, you'd get also
> more numerical stability.
>
> I tried to compile the SVN version a few times the last couple of days, but
> it always broke down: Can someone do a "make clean" in their sandbox and
> repair the error?
>
>  g++ -fPIC -L/home/a907hart/local/lib -L/usr/lib  gdaldem.o commonutils.o
> -L/home/a907hart/src/gdal -lgdal  -larmadillo -L/home/a907hart/local/lib
> -lgeos_c -L/home/a907hart/local/pgsql/lib -lpq -lpthread -lm -lrt -ldl
> -L/home/a907hart/local/lib   -lcurl            -o gdaldem
> gdaldem.o:(.data.rel.ro._ZTV21GDALGeneric3x3Dataset[vtable for
> GDALGeneric3x3Dataset]+0x30): undefined reference to
> `GDALDataset::CloseDependentDatasets()'
> gdaldem.o:(.data.rel.ro._ZTV22GDALColorReliefDataset[vtable for
> GDALColorReliefDataset]+0x30): undefined reference to
> `GDALDataset::CloseDependentDatasets()'
>
>
> Jan
>
> On Thu, Sep 8, 2011 at 5:23 AM, Jan Hartmann <[hidden email]> wrote:
>
> Not sure whether this can be considered a bug, so I give it for what it
> is worth.
>
> I'm doing thin plate spline transformation from one set of projected
> coordinates to another. Both sets have values between -600000 and 600000
> (meters). A typical set of gcps looks like:
>
>  -gcp 62402 -74383 181915 315371 \
>  -gcp -100169 -24213 19685 366661 \
>  -gcp 118323 233822 239994 623190 \
> ...
>
> When transforming the the first two columns of this list with
> gdaltransform -tps, using the same gcp-list, the results should be
> exactly equal to the last two columns:
>
>  (from gdal_tps.cpp:)
>  * The thin plate spline transformer produces exact transformation
>  * at all control points and smoothly varying transformations between
>  * control points with greatest influence from local control points.
>  * It is suitable for for many applications not well modelled by
> polynomial * transformations.
>  *
>
>
> However, they aren't. With a few control points the differences are only
> in millimeters, but with a few hundred gcps the differences become
> meters, and above thousand points the error can get to fifty meters. The
> problem gets worse when some control points are very near to other ones.
> I was happy to discover that dividing all numbers by 1000 solves the
> problem, but there certainly is a numerical instability in the
> algorithm.
>
> Of course, thin plate spline transformations normally use scan
> coordinates as input, and these are almost always small numbers. Even
> so, I can imagine problems with very large rasters and the rubber
> sheeting applications I am working with.
>
> Jan
>
>
>
> _______________________________________________
> gdal-dev mailing list
> [hidden email]
> http://lists.osgeo.org/mailman/listinfo/gdal-dev
>
> _______________________________________________
> gdal-dev mailing list
> [hidden email]
> http://lists.osgeo.org/mailman/listinfo/gdal-dev
>
> _______________________________________________
> gdal-dev mailing list
> [hidden email]
> http://lists.osgeo.org/mailman/listinfo/gdal-dev
>
_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: adding georeferencing information

Jan Hartmann
In reply to this post by Jan Hartmann
Found the answer to this one: use the "listgeo" and "geotifcp" utilities from the Geotiff library (http://trac.osgeo.org/geotiff/). They are not included in the gdal utilities, so you need them install them yourself. There is a debian package for it (geotiff-bin), and I guess there are binary versions for other distributions. Alternatively, you could use fwtools or osgeo4w.

Jan

On 09/05/2011 04:46 PM, Jan Hartmann wrote:
Hi all,

Is there a quick way to add back georeferencing information to a geotif file that has been processed by Imagemagick? I have lots of georeferenced images that I want to sharpen up a bit, but ImageMagick throws the georeference headers away. Is there a way to transport them from the original file to the derived one?

Jan


_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev

_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: adding georeferencing information

Jean-Claude Repetto
Le 22/09/2011 13:37, Jan Hartmann a écrit :
> Found the answer to this one: use the "listgeo" and "geotifcp" utilities
> from the Geotiff library (http://trac.osgeo.org/geotiff/). They are not
> included in the gdal utilities, so you need them install them yourself.
>
> Jan

There is a new GDAL utility program, that has been added recently into
trunk by Even, for that purpose.
It works for ECW files, but I don't know if it works with Tiff files.
Its name is gdal_edit.py :

Usage: gdal_edit [--help-general] [-a_srs srs_def] [-a_ullr ulx uly lrx lry]
                  [-tr xres yres] [-unsetgt] [-a_nodata value]
                  [-mo "META-TAG=VALUE"]*  datasetname

Edit in place various information of an existing GDAL dataset.



Jean-Claude
_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: adding georeferencing information

Jan Hartmann
Thanks Jean-Claude, I'll try it out as soon as I manage to compile GDAL trunk. It's still broken at the moment ...

Jan

On 09/22/2011 02:26 PM, Jean-Claude Repetto wrote:
Le 22/09/2011 13:37, Jan Hartmann a écrit :
Found the answer to this one: use the "listgeo" and "geotifcp" utilities
from the Geotiff library (http://trac.osgeo.org/geotiff/). They are not
included in the gdal utilities, so you need them install them yourself.

Jan

There is a new GDAL utility program, that has been added recently into trunk by Even, for that purpose.
It works for ECW files, but I don't know if it works with Tiff files.
Its name is gdal_edit.py :

Usage: gdal_edit [--help-general] [-a_srs srs_def] [-a_ullr ulx uly lrx lry]
                 [-tr xres yres] [-unsetgt] [-a_nodata value]
                 [-mo "META-TAG=VALUE"]*  datasetname

Edit in place various information of an existing GDAL dataset.



Jean-Claude
_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev

_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: adding georeferencing information

Jan Hartmann
In reply to this post by Jan Hartmann
Hi Andrew, I must have overlooked your mail. Applygeo is not mentioned on the libgeotiff page, but I found the following posting:

http://osgeo-org.1803224.n2.nabble.com/Copy-projection-between-GeoTIFFs-td2065142.html

Thanks,

Jan   

On 09/22/2011 04:24 PM, Andrew Brooks wrote:
On 22 September 2011 12:37, Jan Hartmann <[hidden email]> wrote:
Found the answer to this one: use the "listgeo" and "geotifcp" utilities from the Geotiff library (http://trac.osgeo.org/geotiff/).


Did you not receive the email I sent to the list on 6 September?

I said you should use the applygeo" utility from the libgeotiff distribution, it is better than the listgeo/geotifcp combination.

Andrew



_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: adding georeferencing information

Andrew Brooks
On 22 September 2011 16:06, Jan Hartmann <[hidden email]> wrote:
Hi Andrew, I must have overlooked your mail. Applygeo is not mentioned on the libgeotiff page, but I found the following posting:

http://osgeo-org.1803224.n2.nabble.com/Copy-projection-between-GeoTIFFs-td2065142.html

It is in the bin directory of the libgeotiff distribution:


No idea why it's not documented or mentioned on the libgeotiff web page! It writes the geo information into the TIFF without having to re-write (copy) the TIFF itself.

Andrew


_______________________________________________
gdal-dev mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/gdal-dev
Loading...