Re: [gdal-dev] Gdal and Google's OSS Fuzzing project

[Jesse McGraw]

Nice to see that I'm only 3 weeks behind the curve!

> On May 8, 2017, at 2:58 PM, Kurt Schwehr <schw...@gmail.com> wrote:
> 
> Yup... https://lists.osgeo.org/pipermail/gdal-dev/2017-April/046495.html
> 
> I'd be happy if anyone else wanted to take lead on it.
> 
> I've created https://trac.osgeo.org/gdal/ticket/6883
> 
> Since I'm internal to google, I've been running some fuzzer targets against 
> gdal behind the scenes and used the results to fix a number of bugs.  I've 
> added a number of fuzz targets to 
> https://github.com/schwehr/gdal-autotest2/tree/master/cpp and modified GDAL 
> to make fuzzing more productive... e.g. 
> 
> https://trac.osgeo.org/gdal/changeset/37592/ adds 
> FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION to a driver
> https://trac.osgeo.org/gdal/changeset/37909 example fix
> 
> I have ~50 bugs that I haven't gotten to.
> 
>> On Mon, May 8, 2017 at 11:46 AM, Jesse McGraw <jlmcg...@gmail.com> wrote:
>> I think the gdal suite would be a perfect candidate for this project from 
>> google.  Is anyone interested in trying to integrate gdal into it?
>> 
>> https://opensource.googleblog.com/2017/05/oss-fuzz-five-months-later-and.html?m=1
>> ___
>> gdal-dev mailing list
>> gdal-dev@lists.osgeo.org
>> https://lists.osgeo.org/mailman/listinfo/gdal-dev
> 
> 
> 
> -- 
> --
> http://schwehr.org
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/gdal-dev

[gdal-dev] Gdal and Google's OSS Fuzzing project

[Jesse McGraw]

I think the gdal suite would be a perfect candidate for this project from 
google.  Is anyone interested in trying to integrate gdal into it?

https://opensource.googleblog.com/2017/05/oss-fuzz-five-months-later-and.html?m=1
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/gdal-dev

Re: [gdal-dev] Core dump using gdal_translate to create geopackage

[Jesse McGraw]

Thanks Even.

I've uploaded it to dropbox here: 
https://www.dropbox.com/s/7tcxh9bvk4127c6/Washington_SEC.tif?dl=0




On 02/15/2016 09:53 AM, Even Rouault wrote:

Jesse,

The geotiff would be necessary to replicate.

Even

Le lundi 15 février 2016 15:46:15, Jesse McGraw a écrit :

I was fiddling with Even's geopackage examples
<http://erouault.blogspot.com/2014/12/gdal-geopackage-raster-support.html>
when I came across a core dump issue in gdal_translate

gdal version is compiled from github trunk as of 15 February 2016 under
Ubuntu 15.10 x64

 gdal_translate --version
 GDAL 2.1.0dev, released 2015/99/99



I can share the geotiff in question if need be.  Let me know if I can
provide any other helpful information

Basic command is:

 gdal_translate \
  Washington_SEC.tif \
  test.gpkg \
  -of GPKG \
  -co TILE_FORMAT=PNG \
  -co TILING_SCHEME=GoogleMapsCompatible


Out from Valgrind is:

 jlmcgraw@UbuntuMate1504:~/Documents/myPrograms/gdal_tests$
 CPL_DEBUG=ON \

  > valgrind --tool=memcheck --leak-check=yes --read-var-info=yes \
  > gdal_translate \
  >
  > Washington_SEC.tif \
  > test.gpkg \
  > -of GPKG \
  > -co TILE_FORMAT=PNG \
  > -co TILING_SCHEME=GoogleMapsCompatible

 ==30978== Memcheck, a memory error detector
 ==30978== Copyright (C) 2002-2015, and GNU GPL'd, by Julian Seward
 et al.
 ==30978== Using Valgrind-3.11.0 and LibVEX; rerun with -h for
 copyright info
 ==30978== Command: gdal_translate Washington_SEC.tif test.gpkg -of
 GPKG -co TILE_FORMAT=PNG -co TILING_SCHEME=GoogleMapsCompatible
 ==30978==
 GDAL: GDALOpen(Washington_SEC.tif, this=0xf7e9d80) succeeds as GTiff.
 Input file size is 15538, 11309
 GDAL: GDAL_CACHEMAX = 100 MB
 GDAL: GDALWarpKernel()::GWKGeneralCase()
 Src=0,0,1944x1415 Dst=0,0,2737x1992
 0GDAL: GDALWarpKernel()::GWKGeneralCase()
 Src=1941,0,1945x1415 Dst=2737,0,2737x1992
 .GDAL: GDALWarpKernel()::GWKGeneralCase()
 Src=3883,0,1945x1415 Dst=5474,0,2737x1992
 --30978-- warning: evaluate_Dwarf3_Expr: unhandled DW_OP_ 0xf2
 --30978-- warning: evaluate_Dwarf3_Expr: unhandled DW_OP_ 0xf2
 --30978-- warning: evaluate_Dwarf3_Expr: unhandled DW_OP_ 0xf2
 --30978-- warning: evaluate_Dwarf3_Expr: unhandled DW_OP_ 0xf2
 ==30978== Invalid write of size 8
 ==30978==at 0x4C3259F: memset (in
 /usr/lib/valgrind/vgpreload_memcheck-amd64-linux.so)
 ==30978==by 0x5692CD8: memset (string3.h:90)
 ==30978==by 0x5692CD8: GDALGeoPackageDataset::ReadTile(int, int,
 unsigned char*, bool*) (gdalgeopackagerasterband.cpp:641)
 ==30978==by 0x5692FA2: GDALGeoPackageDataset::ReadTile(int, int)
 (gdalgeopackagerasterband.cpp:489)
 ==30978==by 0x5694B51: GDALGeoPackageRasterBand::IReadBlock(int,
 int, void*) (gdalgeopackagerasterband.cpp:733)
 ==30978==by 0x5553B15: GDALRasterBand::GetLockedBlockRef(int,
 int, int) (gdalrasterband.cpp:1142)
 ==30978==by 0x5572525: GDALRasterBand::IRasterIO(GDALRWFlag,
 int, int, int, int, void*, int, int, GDALDataType, long long, long
 long, GDALRasterIOExtraArg*) (rasterio.cpp:307)
 ==30978==by 0x557312A:
 GDALDataset::BlockBasedRasterIO(GDALRWFlag, int, int, int, int,
 void*, int, int, GDALDataType, int, int*, long long, long long, long
 long, GDALRasterIOExtraArg*) (rasterio.cpp:2938)
 ==30978==by 0x551B09C: GDALDataset::RasterIO(GDALRWFlag, int,
 int, int, int, void*, int, int, GDALDataType, int, int*, long long,
 long long, long long, GDALRasterIOExtraArg*) (gdaldataset.cpp:1991)
 ==30978==by 0x551B11D: GDALDatasetRasterIO (gdaldataset.cpp:2034)
 ==30978==by 0x56384C4: GDALWarpOperation::WarpRegion(int, int,
 int, int, int, int, int, int, int, int, double, double)
 (gdalwarpoperation.cpp:1370)
 ==30978==by 0x56387C6: GDALWarpOperation::ChunkAndWarpImage(int,
 int, int, int) (gdalwarpoperation.cpp:720)
 ==30978==by 0x57A9F5A: GDALGeoPackageDataset::CreateCopy(char
 const*, GDALDataset*, int, char**, int (*)(double, char const*,
 void*), void*) (ogrgeopackagedatasource.cpp:3242)
 ==30978==  Address 0x0 is not stack'd, malloc'd or (recently) free'd
 ==30978==
 ==30978==
 ==30978== Process terminating with default action of signal 11
(SIGSEGV) ==30978==  Access not within mapped region at address 0x0
 ==30978==at 0x4C3259F: memset (in
 /usr/lib/valgrind/vgpreload_memcheck-amd64-linux.so)
 ==30978==by 0x5692CD8: memset (string3.h:90)
 ==30978==by 0x5692CD8: GDALGeoPackageDataset::ReadTile(int, int,
 unsigned char*, bool*) (gdalgeopackagerasterband.cpp:641)
 ==30978==by 0x5692FA2: GDALGeoPackageDataset::ReadTile(int, int)
 (gdalgeopackagerasterband.cpp:489)
 ==30978==by 0x5694B51: GDALGeoPackageRasterBand::

[gdal-dev] Core dump using gdal_translate to create geopackage

[Jesse McGraw]

I was fiddling with Even's geopackage examples 
 
when I came across a core dump issue in gdal_translate


gdal version is compiled from github trunk as of 15 February 2016 under 
Ubuntu 15.10 x64


   gdal_translate --version
   GDAL 2.1.0dev, released 2015/99/99



I can share the geotiff in question if need be.  Let me know if I can 
provide any other helpful information


Basic command is:

   gdal_translate \
Washington_SEC.tif \
test.gpkg \
-of GPKG \
-co TILE_FORMAT=PNG \
-co TILING_SCHEME=GoogleMapsCompatible


Out from Valgrind is:

   jlmcgraw@UbuntuMate1504:~/Documents/myPrograms/gdal_tests$
   CPL_DEBUG=ON \
> valgrind --tool=memcheck --leak-check=yes --read-var-info=yes \
> gdal_translate \
> Washington_SEC.tif \
> test.gpkg \
> -of GPKG \
> -co TILE_FORMAT=PNG \
> -co TILING_SCHEME=GoogleMapsCompatible
   ==30978== Memcheck, a memory error detector
   ==30978== Copyright (C) 2002-2015, and GNU GPL'd, by Julian Seward
   et al.
   ==30978== Using Valgrind-3.11.0 and LibVEX; rerun with -h for
   copyright info
   ==30978== Command: gdal_translate Washington_SEC.tif test.gpkg -of
   GPKG -co TILE_FORMAT=PNG -co TILING_SCHEME=GoogleMapsCompatible
   ==30978==
   GDAL: GDALOpen(Washington_SEC.tif, this=0xf7e9d80) succeeds as GTiff.
   Input file size is 15538, 11309
   GDAL: GDAL_CACHEMAX = 100 MB
   GDAL: GDALWarpKernel()::GWKGeneralCase()
   Src=0,0,1944x1415 Dst=0,0,2737x1992
   0GDAL: GDALWarpKernel()::GWKGeneralCase()
   Src=1941,0,1945x1415 Dst=2737,0,2737x1992
   .GDAL: GDALWarpKernel()::GWKGeneralCase()
   Src=3883,0,1945x1415 Dst=5474,0,2737x1992
   --30978-- warning: evaluate_Dwarf3_Expr: unhandled DW_OP_ 0xf2
   --30978-- warning: evaluate_Dwarf3_Expr: unhandled DW_OP_ 0xf2
   --30978-- warning: evaluate_Dwarf3_Expr: unhandled DW_OP_ 0xf2
   --30978-- warning: evaluate_Dwarf3_Expr: unhandled DW_OP_ 0xf2
   ==30978== Invalid write of size 8
   ==30978==at 0x4C3259F: memset (in
   /usr/lib/valgrind/vgpreload_memcheck-amd64-linux.so)
   ==30978==by 0x5692CD8: memset (string3.h:90)
   ==30978==by 0x5692CD8: GDALGeoPackageDataset::ReadTile(int, int,
   unsigned char*, bool*) (gdalgeopackagerasterband.cpp:641)
   ==30978==by 0x5692FA2: GDALGeoPackageDataset::ReadTile(int, int)
   (gdalgeopackagerasterband.cpp:489)
   ==30978==by 0x5694B51: GDALGeoPackageRasterBand::IReadBlock(int,
   int, void*) (gdalgeopackagerasterband.cpp:733)
   ==30978==by 0x5553B15: GDALRasterBand::GetLockedBlockRef(int,
   int, int) (gdalrasterband.cpp:1142)
   ==30978==by 0x5572525: GDALRasterBand::IRasterIO(GDALRWFlag,
   int, int, int, int, void*, int, int, GDALDataType, long long, long
   long, GDALRasterIOExtraArg*) (rasterio.cpp:307)
   ==30978==by 0x557312A:
   GDALDataset::BlockBasedRasterIO(GDALRWFlag, int, int, int, int,
   void*, int, int, GDALDataType, int, int*, long long, long long, long
   long, GDALRasterIOExtraArg*) (rasterio.cpp:2938)
   ==30978==by 0x551B09C: GDALDataset::RasterIO(GDALRWFlag, int,
   int, int, int, void*, int, int, GDALDataType, int, int*, long long,
   long long, long long, GDALRasterIOExtraArg*) (gdaldataset.cpp:1991)
   ==30978==by 0x551B11D: GDALDatasetRasterIO (gdaldataset.cpp:2034)
   ==30978==by 0x56384C4: GDALWarpOperation::WarpRegion(int, int,
   int, int, int, int, int, int, int, int, double, double)
   (gdalwarpoperation.cpp:1370)
   ==30978==by 0x56387C6: GDALWarpOperation::ChunkAndWarpImage(int,
   int, int, int) (gdalwarpoperation.cpp:720)
   ==30978==by 0x57A9F5A: GDALGeoPackageDataset::CreateCopy(char
   const*, GDALDataset*, int, char**, int (*)(double, char const*,
   void*), void*) (ogrgeopackagedatasource.cpp:3242)
   ==30978==  Address 0x0 is not stack'd, malloc'd or (recently) free'd
   ==30978==
   ==30978==
   ==30978== Process terminating with default action of signal 11 (SIGSEGV)
   ==30978==  Access not within mapped region at address 0x0
   ==30978==at 0x4C3259F: memset (in
   /usr/lib/valgrind/vgpreload_memcheck-amd64-linux.so)
   ==30978==by 0x5692CD8: memset (string3.h:90)
   ==30978==by 0x5692CD8: GDALGeoPackageDataset::ReadTile(int, int,
   unsigned char*, bool*) (gdalgeopackagerasterband.cpp:641)
   ==30978==by 0x5692FA2: GDALGeoPackageDataset::ReadTile(int, int)
   (gdalgeopackagerasterband.cpp:489)
   ==30978==by 0x5694B51: GDALGeoPackageRasterBand::IReadBlock(int,
   int, void*) (gdalgeopackagerasterband.cpp:733)
   ==30978==by 0x5553B15: GDALRasterBand::GetLockedBlockRef(int,
   int, int) (gdalrasterband.cpp:1142)
   ==30978==by 0x5572525: GDALRasterBand::IRasterIO(GDALRWFlag,
   int, int, int, int, void*, int, int, GDALDataType, long long, long
   long, GDALRasterIOExtraArg*) (rasterio.cpp:307)
   ==30978==by 0x557312A:
   GDALDataset::BlockBasedRasterIO(GDALRWFlag, int, int, int, int,
   void*, int, 

Re: [gdal-dev] gdalwarp -crop_to_cutline produces non-square pixels for EPSG:3857 (GDAL 1.10.1, released 2013/08/26)

[Jesse McGraw]

Thanks Even.

I did try re-warping the cropped output to 3857 again and it does produce
square pixels but at the expense of an additional warping operation

I'm just curious why the warping that is part of the -crop_to_cutline
process is changing the pixel ratio while a regular, non-cropping warping
makes them square

On Thu, Feb 5, 2015 at 11:36 AM, Even Rouault even.roua...@spatialys.com
wrote:

 Le jeudi 05 février 2015 17:23:06, Jesse McGraw a écrit :
  I'm not sure whether this is a real issue or not but I thought I'd bring
 it
  up, at the very least I'll learn something
 
  When using gdalwarp -cutline shapefile -crop_to_cutline on an input
  raster that is in EPSG:3857 with square-pixels the output raster, while
  still EPSG:3857, now has non-square pixels.
 
  They're aren't terribly non-square but aren't they supposed to be
  completely square for EPSG:3857?

 Perhaps... Actually the effet of -crop_to_cutline is exactly the same as
 manually passing the target extents with -te with the bounding box of the
 cutline.
 It is not possible to both preserve the extent extents and pixel size at
 the
 same time, due to width and height being integer values.
 So, in order to preserve pixel square, the extent should be modified a
 little
 bit. What is more appropriate is a matter of point of view I think.
 You can always re-run gdalwarp tmp.tif out.tif where tmp.tif is the
 output
 of first gdalwarp with -crop_to_cutline. And you should get square pixels I
 believe.

 
  (FWIW, I see that there are tickets opened that reference similar issues
  but they reference the output being shifted or the origin changing, not
 the
  pixel shape changing)
 
  For example:
  #Warp our original .tif to EPSG:3857
  $gdalwarp \
  -t_srs EPSG:3857 \
  -dstalpha \
  -co TILED=YES \
  ENR_L33.tif \
  ./2.tif
 
  #See that the output pixels are square
  $gdalinfo 2.tif
  Origin = (-8577554.996301921084523,5421778.172851986251771)
  Pixel Size = (43.677179501975118,-43.677179501975118)
 
 
  #Now crop the image to a cutline
  $gdalwarp \
  -crop_to_cutline \
  -dstalpha \
  -cutline ./ENR_L33.shp \
  -cblend 10 \
  -co TILED=YES \
  ./2.tif \
  ./3.tif
 
  #See that output pixels are not square
  $ gdalinfo 3.tif
  Origin = (-8480047.445924906060100,5366376.137789577245712)
  Pixel Size = (43.678439570399853,-43.675457269061347)
 
 
  I ran some more tests and without the -crop_to_cutline option the output
  pixels remain square
 
  Thanks,
Jesse

 --
 Spatialys - Geospatial professional services
 http://www.spatialys.com

___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/gdal-dev

[gdal-dev] gdalwarp -crop_to_cutline produces non-square pixels for EPSG:3857 (GDAL 1.10.1, released 2013/08/26)

[Jesse McGraw]

I'm not sure whether this is a real issue or not but I thought I'd bring it
up, at the very least I'll learn something

When using gdalwarp -cutline shapefile -crop_to_cutline on an input
raster that is in EPSG:3857 with square-pixels the output raster, while
still EPSG:3857, now has non-square pixels.

They're aren't terribly non-square but aren't they supposed to be
completely square for EPSG:3857?

(FWIW, I see that there are tickets opened that reference similar issues
but they reference the output being shifted or the origin changing, not the
pixel shape changing)

For example:
#Warp our original .tif to EPSG:3857
$gdalwarp \
-t_srs EPSG:3857 \
-dstalpha \
-co TILED=YES \
ENR_L33.tif \
./2.tif

#See that the output pixels are square
$gdalinfo 2.tif
Origin = (-8577554.996301921084523,5421778.172851986251771)
Pixel Size = (43.677179501975118,-43.677179501975118)


#Now crop the image to a cutline
$gdalwarp \
-crop_to_cutline \
-dstalpha \
-cutline ./ENR_L33.shp \
-cblend 10 \
-co TILED=YES \
./2.tif \
./3.tif

#See that output pixels are not square
$ gdalinfo 3.tif
Origin = (-8480047.445924906060100,5366376.137789577245712)
Pixel Size = (43.678439570399853,-43.675457269061347)


I ran some more tests and without the -crop_to_cutline option the output
pixels remain square

Thanks,
  Jesse
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/gdal-dev

Re: [gdal-dev] gdalwarp -crop_to_cutline produces non-square pixels for EPSG:3857 (GDAL 1.10.1, released 2013/08/26)

[Jesse McGraw]

Thank you, Even.  For now I'll proceed with the additional warping.

Just for my education: are square pixels a requirement for EPSG:3857?

On Thu, Feb 5, 2015 at 11:50 AM, Even Rouault even.roua...@spatialys.com
wrote:

 Le jeudi 05 février 2015 17:41:22, Jesse McGraw a écrit :
  Thanks Even.
 
  I did try re-warping the cropped output to 3857 again and it does produce
  square pixels but at the expense of an additional warping operation
 
  I'm just curious why the warping that is part of the -crop_to_cutline
  process is changing the pixel ratio while a regular, non-cropping warping
  makes them square

 I think I explained it in my previous answer. This is because
 -crop_to_cutline
 uses -te internally, which in the general case doesn't allow pixel squares.
 And non-cropping warping doesn't necessarily preserve extent of the output
 to
 be the same as extent of the input (even in non-reprojection cases). You
 can
 check that actually by comparing gdalinfo on the temporary and final
 datasets.
 There's no fundamental reason why -crop_to_cutline couldn't produce exact
 square pixels, if people are annoyed by that behaviour. Just someone has
 to do
 it...

 
  On Thu, Feb 5, 2015 at 11:36 AM, Even Rouault 
 even.roua...@spatialys.com
 
  wrote:
   Le jeudi 05 février 2015 17:23:06, Jesse McGraw a écrit :
I'm not sure whether this is a real issue or not but I thought I'd
bring
  
   it
  
up, at the very least I'll learn something
   
When using gdalwarp -cutline shapefile -crop_to_cutline on an
 input
raster that is in EPSG:3857 with square-pixels the output raster,
 while
still EPSG:3857, now has non-square pixels.
   
They're aren't terribly non-square but aren't they supposed to be
completely square for EPSG:3857?
  
   Perhaps... Actually the effet of -crop_to_cutline is exactly the same
 as
   manually passing the target extents with -te with the bounding box of
 the
   cutline.
   It is not possible to both preserve the extent extents and pixel size
 at
   the
   same time, due to width and height being integer values.
   So, in order to preserve pixel square, the extent should be modified a
   little
   bit. What is more appropriate is a matter of point of view I think.
   You can always re-run gdalwarp tmp.tif out.tif where tmp.tif is the
   output
   of first gdalwarp with -crop_to_cutline. And you should get square
 pixels
   I believe.
  
(FWIW, I see that there are tickets opened that reference similar
issues but they reference the output being shifted or the origin
changing, not
  
   the
  
pixel shape changing)
   
For example:
#Warp our original .tif to EPSG:3857
$gdalwarp \
   
-t_srs EPSG:3857 \
-dstalpha \
-co TILED=YES \
ENR_L33.tif \
./2.tif
   
#See that the output pixels are square
$gdalinfo 2.tif
Origin = (-8577554.996301921084523,5421778.172851986251771)
Pixel Size = (43.677179501975118,-43.677179501975118)
   
   
#Now crop the image to a cutline
$gdalwarp \
   
-crop_to_cutline \
-dstalpha \
-cutline ./ENR_L33.shp \
-cblend 10 \
-co TILED=YES \
./2.tif \
./3.tif
   
#See that output pixels are not square
$ gdalinfo 3.tif
Origin = (-8480047.445924906060100,5366376.137789577245712)
Pixel Size = (43.678439570399853,-43.675457269061347)
   
   
I ran some more tests and without the -crop_to_cutline option the
output pixels remain square
   
Thanks,
   
  Jesse
  
   --
   Spatialys - Geospatial professional services
   http://www.spatialys.com

 --
 Spatialys - Geospatial professional services
 http://www.spatialys.com

___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/gdal-dev

Re: [gdal-dev] gdalwarp in MapTiler Pro Command Line

[Jesse McGraw]

I have done some work on this using only gdal and other open-source tools,
please see my github repository for ideas:
https://github.com/jlmcgraw?tab=repositories
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/gdal-dev

[gdal-dev] Forcing gdal translate to only rotate/scale

[Jesse McGraw]

Is it possible to force gdal_translate to only rotate/scale an image
when given a list of more than two GCPs?

I'm manually georeferencing maps which are all North-up and I'd like
to supply as many GCPs as possible to average out my errors but if I
provide more than two it will potentially skew the image.

Thanks!
-Jesse
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/gdal-dev

Re: [gdal-dev] gdaltranslate and gdalwarp rotating and resizing image without being asked to, why?

[Jesse McGraw]

Even,
  Thank you for the very detailed responses.  Time to go do some more
googling on GDAL perl bindings (or dig into python)



On Tue, Jul 1, 2014 at 4:37 PM, Even Rouault even.roua...@mines-paris.org
wrote:

 Le mardi 01 juillet 2014 15:50:18, Jesse McGraw a écrit :
  On Mon, Jun 30, 2014 at 3:07 PM, Even Rouault
  even.roua...@mines-paris.org
 
  wrote:
   That might the problem. If you create the GCPs on the PDF and not on
 the
   PNG,
   and that the PDF has rotation, pdftoppm I think will apply the
 rotation,
   hence
   the axis swapping.
 
  Remember, I'm using the PNG I've created for all gdal operations, gdal
  never interacts with the source PDF
 
  I'm probably fundamentally misunderstanding some aspect of the command
  lines I'm using in gdal*.  My end goal is to get the necessary affine
  transformation information to display the PNG I create in step 3a
 properly
  in a separate Android moving-map application (eg. your current position
  displayed on the diagram as you move around the airport)

 Jesse,

 hum ok, with the images I think I now understand the situation (I've only
 looked at the landscape diagram) and what you want.

 You probably don't want to use gdalwarp, which by default will result in a
 north-up image. You likely want to compute the geotransform affine matrix
 for
 your GCPs and set it, since indeed your GCPs define a 1st-order
 transformation.

 The GDAL C API has GDALGCPsToGeoTransform() which does exactly this. I see
 it
 is also available in the Python API :  gdal.GCPsToGeoTransform().
 See http://svn.osgeo.org/gdal/trunk/autotest/gcore/gcps2geotransform.py
 for
 how to use it.

 You could then just do a CreateCopy() of the PNG to a VRT, and set the
 returned geotransform by the previous step with SetGeoTransform().

 Even


 --
 Geospatial professional services
 http://even.rouault.free.fr/services.html

___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/gdal-dev

Re: [gdal-dev] gdaltranslate and gdalwarp rotating and resizing image without being asked to, why?

[Jesse McGraw]

On Mon, Jun 30, 2014 at 3:07 PM, Even Rouault even.roua...@mines-paris.org
wrote:



 That might the problem. If you create the GCPs on the PDF and not on the
 PNG,
 and that the PDF has rotation, pdftoppm I think will apply the rotation,
 hence
 the axis swapping.


Remember, I'm using the PNG I've created for all gdal operations, gdal
never interacts with the source PDF

I'm probably fundamentally misunderstanding some aspect of the command
lines I'm using in gdal*.  My end goal is to get the necessary affine
transformation information to display the PNG I create in step 3a properly
in a separate Android moving-map application (eg. your current position
displayed on the diagram as you move around the airport)

Well I'm really confused by what you are trying to do.

 For a better understanding, you should perhaps provide the PDF and all
 exact
 command lines you are typing.


Here's a step-by-step with links to relevant files for a landscape diagram
(SFO)

1. Get source PDFs from the Federal Aviation Administrator (FAA)

https://drive.google.com/file/d/0B4EbUsrclSoDeXdmM3gyV29WN0U/edit?usp=sharing

2. Run my utility.

3a. It converts the PDF to a PNG,

https://drive.google.com/file/d/0B4EbUsrclSoDMTlDZzQzeDJWU1U/edit?usp=sharing

3b. create GCPs for GDAL based on the coordinates within the PNG (hard to
see but green circles mark the GCPs chosen to correspond to the lat/lon
shown on diagram)

https://drive.google.com/file/d/0B4EbUsrclSoDWUtqNVgtWFUwa0k/edit?usp=sharing

3c.  gdal_translate the PNG file based on the GCPs we found to an
intermediate VRT
gdal_translate -q -of VRT -strict -a_srs EPSG:4326  -gcp 1017.14863049096
1800.51717680381 -122.4 37.61667 -gcp 337.854728682168 1260.25
-122.3833 37.6 -gcp 1017.14863049095 1260.25
-122.3833 37.61667 -gcp 1017.14863049096 720.1667
-122.3667 37.61667 -gcp 337.854728682172 720.1667
-122.3667 37.6 -gcp 337.854728682172 1800.33895140302
-122.4 37.6 './dtpp/00375AD.png'
'./dtpp/CA-SFO-00375AD-PDF-AIRPORT-DIAGRAM.vrt'


https://drive.google.com/file/d/0B4EbUsrclSoDX2hzTTBsY1c5dXM/edit?usp=sharing

3d. gdalwarp the intermediate VRT to the final VRT
gdalwarp -q -of VRT -t_srs EPSG:4326 -order 1 -overwrite
''./dtpp/CA-SFO-00375AD-PDF-AIRPORT-DIAGRAM.vrt''
'./dtpp/warpedCA-SFO-00375AD.PDF-AIRPORT DIAGRAM.vrt'


https://drive.google.com/file/d/0B4EbUsrclSoDTmNDanFmbkVvVm8/edit?usp=sharing

Here's the process for a portrait one

1. Get source PDFs from the Federal Aviation Administrator (FAA)

https://drive.google.com/file/d/0B4EbUsrclSoDMjNPOVBPT01JS1U/edit?usp=sharing

2. Run my utility.

3a. It converts the PDF to a PNG,

https://drive.google.com/file/d/0B4EbUsrclSoDTUZNY2RKaVVHMDQ/edit?usp=sharing

3b. create GCPs for GDAL based on the coordinates within the PNG (hard to
see but green circles mark the GCPs chosen to correspond to the lat/lon
shown on diagram)

https://drive.google.com/file/d/0B4EbUsrclSoDbHNyTnhGczh3RDA/edit?usp=sharing

3c.  gdal_translate the PNG file based on the GCPs we found to an
intermediate VRT
 gdal_translate -q -of VRT -strict -a_srs EPSG:4326  -gcp 310.408626924008
1496.72390271202 -77.3 37.5 -gcp 976.562655349203
660.697730425566 -77.31667 37.51667 -gcp 310.523817613394
660.744792369168 -77.3 37.51667 -gcp 976.580225174465
1496.63921742496 -77.31667 37.5 './dtpp/00347AD.png'
'./dtpp/VA-RIC-00347AD-PDF-AIRPORT-DIAGRAM.vrt'


https://drive.google.com/file/d/0B4EbUsrclSoDRzdVZDQzY2lZT1U/edit?usp=sharing

3d. gdalwarp the intermediate VRT to the final VRT
gdalwarp -q -of VRT -t_srs EPSG:4326 -order 1 -overwrite
''./dtpp/VA-RIC-00347AD-PDF-AIRPORT-DIAGRAM.vrt''
'./dtpp/warpedVA-RIC-00347AD.PDF-AIRPORT DIAGRAM.vrt'


https://drive.google.com/file/d/0B4EbUsrclSoDRGptZFhhaVdaT1E/edit?usp=sharing
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/gdal-dev

Re: [gdal-dev] gdaltranslate and gdalwarp rotating and resizing image without being asked to, why?

[Jesse McGraw]

Even,
  Thank you for the response!

I'm actually rasterizing the PDF to PNG via pdftoppm before doing anything
with gdal*, that's why the rotating and resizing surprises me.

When do the pixelSkew parameters get used vs. swapping the axes?  The
swapping works in this case is because the rotation is 90 degrees but what
would happen with an arbitrary rotation (eg 37 degrees).

The reason I'm curious about all this is that I'm wanting to use only the
initially generated PNG for my final georeferenced image and not have to
generate an additional one that is rotated by gdal with the skew parameters
set to 0 and also resized.

If I can provide any more detail on my process or what I'm trying to
accomplish please let me know.

Here's the general flow for Airport Diagrams:

Source PDF - create PNG, create list of GCPs from the PDF mapped to that
PNG - gdal_translate PNG to intermediate VRT - gdal_warp intermediate VRT
to final VRT

The final VRT displays properly in QGIS so the process is working correctly
but I was hoping to be able to extract from it all of the parameters of the
affine transformation necessary to properly display the PNG as a
geoferenced image (origin, extents, pixel sizes and pixel skew) in a
separate Android application.  (This bitmap is to be used in a separate
moving map application to display an airplane's real world position on the
image)

-Jesse


On Sat, Jun 28, 2014 at 12:01 PM, Even Rouault even.roua...@mines-paris.org
 wrote:

 Jesse,

 I'm not completely sure to really understand the axis swapping issue
 without
 reproducing what you do.
 The PDF driver should take into account the PDF rotation parameter when
 rasterizing it (so that it outputs something that is visually identical
 with
 what other PDF viewers do).

 To answer your questions:

 1) When not providing an existing target image, gdalwarp will indeed
 rectify
 images that have an arbitrary transform with skew/rotation terms, or GCPs,
 to
 an affine transform with squared pixels, which is generally the desired
 behaviour. If you prepare yourself the target dataset with the geotransform
 you whish, you should be able run gdalwarp on it. You can alternatively try
 the -ts option of gdalwarp to force the target dimensions.

 2) The resizing is a consequence of the output geotransform being different
 that the source geotransform.

 3) You could likely merge your step 1 and 2, i.e. directly generate a .png
 with the SRS and GCP set. However I'd recommand generating a tiled TIFF
 instead to make gdalwarp faster.

 Best regards,

 Even

 Le vendredi 27 juin 2014 19:49:35, Jesse McGraw a écrit :
  Gdal gurus,
 
  I'm writing a utility to automatically georeference Instrument Approach
  Procedures and Airport Diagrams from the FAA (github.com/jlmcgraw, if
  you're curious) and I'm hoping for some insight from the experts here as
 to
  why I'm seeing the following behavior with Airport Diagrams
 
  Most of the Airport Diagrams diagrams are portrait orientation (True
 North
  up), however some of them are actually laid out in landscape orientation
  (True North to the left or right) though the source PDF is still in
  portrait
 
  If I georeference these landscape oriented using GCPs then gdalwarp will
  swap X and Y axes and also resize the resulting bitmap (it also resizes
  portrait oriented ones too, though without swapping axes)
 
  My questions:
  1) Why is gdalwarp swapping the axes rather than using the X and Y
  pixel-skew parameters of the affine transformation (seems to be always
 0)?
  2) Why is it re-sizing the bitmap?
  3) Do I have to do both of these steps separately (translate and then
 warp)
  or can I simplify somehow?
 
  Steps in the process:
 
  1) Generate a list of Ground Control Points from the PDF, rasterize the
 PDF
  to a 300dpi PNG
 
  2) gdal_translate .PNG to an intermediate .VRT
  (gdal_translate -q -of VRT -strict -a_srs EPSG:4326 $gcpstring
  '$targetpng'  '$targetvrt')
 
  3) gdalwarp intermediate .VRT to get the final georeferenced output .VRT
   (gdalwarp -q -of VRT -t_srs EPSG:4326 -order 1 -overwrite
 ''$targetvrt''
  '$targetvrt2')
 
 
 
 
  *example from a landscape diagram*
 
  (snippet of step 2 output from gdal_translate'd .vrt)
  CA-SFO-00375AD-PDF-AIRPORT-DIAGRAM.vrt:
  *VRTDataset rasterXSize=1613 rasterYSize=2475*
 
  (snippet of step 3 output from gdalwarp'd .vrt)
  warpedCA-SFO-00375AD.PDF-AIRPORT\ DIAGRAM.vrt:
  *VRTDataset rasterXSize=2623 rasterYSize=1359
  subClass=VRTWarpedDataset*
  *GeoTransform -1.2242081739746585e+02,  2.9114844968978991e-05,
  0.e+00,  3.7641622693302779e+01,  0.e+00,
  -2.9114844968978991e-05/GeoTransform*
 
 
  *example from a portrait diagram*
(snippet of step 2 output from gdal_translate'd .vrt)
  *VA-RIC-00347AD-PDF-AIRPORT-*
  *DIAGRAM.vrt*
  *VRTDataset rasterXSize=1613 rasterYSize=2475*
 
  (snippet of step 3 output from gdalwarp'd .vrt)
  *warpedVA-RIC-00347AD.PDF-*
  *AIRPORT\ DIAGRAM.vrt

[gdal-dev] gdaltranslate and gdalwarp rotating and resizing image without being asked to, why?

[Jesse McGraw]

Gdal gurus,

I'm writing a utility to automatically georeference Instrument Approach
Procedures and Airport Diagrams from the FAA (github.com/jlmcgraw, if
you're curious) and I'm hoping for some insight from the experts here as to
why I'm seeing the following behavior with Airport Diagrams

Most of the Airport Diagrams diagrams are portrait orientation (True North
up), however some of them are actually laid out in landscape orientation
(True North to the left or right) though the source PDF is still in portrait

If I georeference these landscape oriented using GCPs then gdalwarp will
swap X and Y axes and also resize the resulting bitmap (it also resizes
portrait oriented ones too, though without swapping axes)

My questions:
1) Why is gdalwarp swapping the axes rather than using the X and Y
pixel-skew parameters of the affine transformation (seems to be always 0)?
2) Why is it re-sizing the bitmap?
3) Do I have to do both of these steps separately (translate and then warp)
or can I simplify somehow?

Steps in the process:

1) Generate a list of Ground Control Points from the PDF, rasterize the PDF
to a 300dpi PNG

2) gdal_translate .PNG to an intermediate .VRT
(gdal_translate -q -of VRT -strict -a_srs EPSG:4326 $gcpstring
'$targetpng'  '$targetvrt')

3) gdalwarp intermediate .VRT to get the final georeferenced output .VRT
 (gdalwarp -q -of VRT -t_srs EPSG:4326 -order 1 -overwrite ''$targetvrt''
'$targetvrt2')




*example from a landscape diagram*

(snippet of step 2 output from gdal_translate'd .vrt)
CA-SFO-00375AD-PDF-AIRPORT-DIAGRAM.vrt:
*VRTDataset rasterXSize=1613 rasterYSize=2475*

(snippet of step 3 output from gdalwarp'd .vrt)
warpedCA-SFO-00375AD.PDF-AIRPORT\ DIAGRAM.vrt:
*VRTDataset rasterXSize=2623 rasterYSize=1359
subClass=VRTWarpedDataset*
*GeoTransform -1.2242081739746585e+02,  2.9114844968978991e-05,
0.e+00,  3.7641622693302779e+01,  0.e+00,
-2.9114844968978991e-05/GeoTransform*


*example from a portrait diagram*
  (snippet of step 2 output from gdal_translate'd .vrt)
*VA-RIC-00347AD-PDF-AIRPORT-*
*DIAGRAM.vrt*
*VRTDataset rasterXSize=1613 rasterYSize=2475*

(snippet of step 3 output from gdalwarp'd .vrt)
*warpedVA-RIC-00347AD.PDF-*
*AIRPORT\ DIAGRAM.vrt*



*VRTDataset rasterXSize=1870 rasterYSize=2287
subClass=VRTWarpedDatasetGeoTransform -7.7341103107298423e+01,
2.1580033523515700e-05,  0.e+00,  3.7529840835409082e+01,
0.e+00, -2.1580033523515700e-05/GeoTransform*
Thanks!
-Jesse McGraw
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/gdal-dev

[gdal-dev] ogr2ogr v1.10.1 and AIXM (issues with curves?)

[Jesse McGraw]

GDAL gurus,
 I'm trying to use the AIXM files from the FAA (downloadable here:
https://nfdc.faa.gov/webContent/56DaySub/2014-05-29/aixm5.0.zip ) in QGIS
and I'm noticing that there are some that with no warnings or complaints
just don't produce a polygon. (ogrinfo says POLYGON EMPTY)

Two examples are P-56B DISTRICT OF COLUMBIA.xml and POINSETT MOA, SC.xml

Is there some command line incantation that I'm missing?  Or perhaps this
is a feature of GML/AIXM that isn't supported?  I've read that there is
little or no GDAL support for curves/arcs though I'm unclear as to the
current status of that support.

Here's a snippet from the P-56B file:
horizontalProjection
 Surface gml:id=Surface1 srsDimension=2
srsName=URN:OGC:DEF:CRS:OGC:1.3:CRS84
  patches
   PolygonPatch
exterior
 Ring
  curveMember
   Curve ns:id=Curve1
segments
 CircleByCenterPoint numArc=1
  pointProperty
   Point ns:id=Point_1
pos-77.066944 38.921389/pos
/Point
   /pointProperty
  radius uom=NM0.4/radius
 /CircleByCenterPoint
/segments
   /Curve
  /curveMember
 /Ring
/exterior
   /PolygonPatch
  /patches
 /Surface
/horizontalProjection

If I should be posting this somewhere else or filing a bug report/feature
request please let me know

Thanks!
-Jesse McGraw
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/gdal-dev

Re: [gdal-dev] ogr2ogr v1.10.1 and AIXM (issues with curves?)

[Jesse McGraw]

Well, that explains that one.  :)

I'm sure it must be the same story for the other files, various unsupported
GML features.  Is there any way to make ogr2ogr warn when it encounters
something it doesn't support versus silently creating blank geometry?

Out of curiosity, is more robust GML/AIXM support something that's of
general interest?  Is it widely used?




On Fri, May 16, 2014 at 2:44 PM, Even Rouault
even.roua...@mines-paris.orgwrote:

 Le vendredi 16 mai 2014 20:18:22, Jesse McGraw a écrit :
  GDAL gurus,
   I'm trying to use the AIXM files from the FAA (downloadable here:
  https://nfdc.faa.gov/webContent/56DaySub/2014-05-29/aixm5.0.zip ) in
 QGIS
  and I'm noticing that there are some that with no warnings or complaints
  just don't produce a polygon. (ogrinfo says POLYGON EMPTY)
 
  Two examples are P-56B DISTRICT OF COLUMBIA.xml and POINSETT MOA,
  SC.xml
 
  Is there some command line incantation that I'm missing? Or perhaps this
  is a feature of GML/AIXM that isn't supported?  I've read that there is
  little or no GDAL support for curves/arcs though I'm unclear as to the
  current status of that support.

 CircleByCenterPoint is not currently supported by the GML decoder.

 
  Here's a snippet from the P-56B file:
  horizontalProjection
   Surface gml:id=Surface1 srsDimension=2
  srsName=URN:OGC:DEF:CRS:OGC:1.3:CRS84
patches
 PolygonPatch
  exterior
   Ring
curveMember
 Curve ns:id=Curve1
  segments
   CircleByCenterPoint numArc=1
pointProperty
 Point ns:id=Point_1
  pos-77.066944 38.921389/pos
  /Point
 /pointProperty
radius uom=NM0.4/radius
   /CircleByCenterPoint
  /segments
 /Curve
/curveMember
   /Ring
  /exterior
 /PolygonPatch
/patches
   /Surface
  /horizontalProjection
 
  If I should be posting this somewhere else or filing a bug report/feature
  request please let me know

 That could be worth one.

 
  Thanks!
  -Jesse McGraw

 --
 Geospatial professional services
 http://even.rouault.free.fr/services.html

___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/gdal-dev