Re: [gdal-dev] layer.GetSpatialRef() fails on linux for shapefiles
Hi Javier,
I've come across some other bugs that are specific to Conda running in
Pycharm only this morning. Up to and including getting different
projection strings back for the exact same file, depending on whether
the `GetProjectionRef` is called via Pycharm or the terminal (using the
exact same conda environment).
`osr.SpatialReference().ImportFromEPSG(4326)` returns a RuntimeError in
Pycharm but works in the terminal:
PROJ: proj_create_from_database: Open of
/some/path/anaconda3/envs/abc/share/proj failed
So there's something weird going on with environment's that's Posix
specific.
Thanks,
Jonathan
On 2023-10-12 09:10, Javier Jimenez Shaw wrote:
Are you running pycharm under the conda environment? If you lose any
setting while combining both it can be a problem.
Have you tried setting PROJ_DEBUG=3 ?
It will show some traces from PROJ, like the location of proj.db used.
The point is if you can see the std out.
On Wed, 11 Oct 2023, 15:18 Jonathan Moules via gdal-dev,
wrote:
Hi List,
So, after more investigation:
* Using Python Anaconda on either Mac or Linux for GDAL (3.7.1),
`GetSpatialRef` triggers a RunTime Error for all shapefiles (but only
shapefiles).
* This happens on Ubuntu (two machines), and a Mac, but only under
PyCharm.
* Using the exact same `conda` environment as triggers the above:
running it directly from the terminal works fine.
So there's something about the combination of conda and Pycharm that
breaks this aspect of GDAL shapefile handling on *nix and Mac.
Any thoughts? Also, who do I actually report this bug to? Is it GDAL,
Conda, Pycharm, something else?
Cheers,
Jonathan
On 2023-09-28 14:58, Jonathan Moules via gdal-dev wrote:
> Well, it seems that PROJ_DATA isn't set in their environment.
But it's
> not set in mine either (`print(os.environ['PROJ_DATA']`)! So no
idea
> why mine works just fine without it (Windows 11 thing?).
>
> Creating a PROJ_DATA env var didn't fix anything. Even adding it
> explicitly in Python.
>
> Their log file does have this in at a WARNING level:
>
> `PROJ: proj_create_from_database: Open of
> /home/user/anaconda3/envs/env1/share/proj failed`
>
> That path has: `drwxrwxr-x` permissions.
>
> To answer Even's earlier question:
>
> `ogrinfo /path/to/shape.shp` works fine on their system.
>
>
>
> On 2023-09-28 12:37, Rahkonen Jukka wrote:
>> Hi,
>>
>> Then they should add that environment if they do not know that
they
>> do not belong to "most users"
>> https://proj.org/en/9.3/usage/environmentvars.html
>>
>> -Jukka Rahkonen-
>>
>> -Alkuperäinen viesti-
>> Lähettäjä: gdal-dev Puolesta
>> Jonathan Moules via gdal-dev
>> Lähetetty: torstai 28. syyskuuta 2023 14.10
>> Vastaanottaja: gdal-dev@lists.osgeo.org; Even Rouault
>>
>> Aihe: Re: [gdal-dev] layer.GetSpatialRef() fails on linux for
shapefiles
>>
>> Hi Even,
>>
>> The colleague doesn't have either a PROJ_LIB or a PROJ_DATA
>> environment variable.
>>
>> I asked another colleague to try it; they're on Ubuntu 20.04,
and it
>> worked for them. I believe using the same setup instructions.
>>
>> Cheers,
>>
>> Jonathan
>>
>> On 2023-09-24 22:37, Jonathan Moules via gdal-dev wrote:
>>> Thanks Even. I don't have access to the machine either as the
>>> colleague is moving to another project. I'll have to see if it
fails
>>> for another *nix user.
>>>
>>> On 2023-09-24 22:35, Even Rouault wrote:
>>>> Le 24/09/2023 à 23:22, Jonathan Moules via gdal-dev a écrit :
>>>>>> Also check if the environment isn't messed up regarding
PROJ and
>>>>> the PROJ_LIB/PROJ_DATA environment variable
>>>>>
>>>>> Thanks Even; sorry, what does this line mean? I'm guessing
you're
>>>>> referring to:
>>>>> https://proj.org/en/9.3/usage/environmentvars.html - what
would a
>>>>> "messed up" one look like?
>>>>>
>>>> Hard to say without access to the machine. Perhaps just try to
>>>> recreate a new Conda env from scratch
>>>>
>>>>
>>> ___
>>> gdal-dev mailing list
>>> gdal-dev@lists.osgeo.org
>>> https://list/
>>> s.osgeo.org
<http://s.osg
Re: [gdal-dev] layer.GetSpatialRef() fails on linux for shapefiles
Hi List, So, after more investigation: * Using Python Anaconda on either Mac or Linux for GDAL (3.7.1), `GetSpatialRef` triggers a RunTime Error for all shapefiles (but only shapefiles). * This happens on Ubuntu (two machines), and a Mac, but only under PyCharm. * Using the exact same `conda` environment as triggers the above: running it directly from the terminal works fine. So there's something about the combination of conda and Pycharm that breaks this aspect of GDAL shapefile handling on *nix and Mac. Any thoughts? Also, who do I actually report this bug to? Is it GDAL, Conda, Pycharm, something else? Cheers, Jonathan On 2023-09-28 14:58, Jonathan Moules via gdal-dev wrote: Well, it seems that PROJ_DATA isn't set in their environment. But it's not set in mine either (`print(os.environ['PROJ_DATA']`)! So no idea why mine works just fine without it (Windows 11 thing?). Creating a PROJ_DATA env var didn't fix anything. Even adding it explicitly in Python. Their log file does have this in at a WARNING level: `PROJ: proj_create_from_database: Open of /home/user/anaconda3/envs/env1/share/proj failed` That path has: `drwxrwxr-x` permissions. To answer Even's earlier question: `ogrinfo /path/to/shape.shp` works fine on their system. On 2023-09-28 12:37, Rahkonen Jukka wrote: Hi, Then they should add that environment if they do not know that they do not belong to "most users" https://proj.org/en/9.3/usage/environmentvars.html -Jukka Rahkonen- -Alkuperäinen viesti- Lähettäjä: gdal-dev Puolesta Jonathan Moules via gdal-dev Lähetetty: torstai 28. syyskuuta 2023 14.10 Vastaanottaja: gdal-dev@lists.osgeo.org; Even Rouault Aihe: Re: [gdal-dev] layer.GetSpatialRef() fails on linux for shapefiles Hi Even, The colleague doesn't have either a PROJ_LIB or a PROJ_DATA environment variable. I asked another colleague to try it; they're on Ubuntu 20.04, and it worked for them. I believe using the same setup instructions. Cheers, Jonathan On 2023-09-24 22:37, Jonathan Moules via gdal-dev wrote: Thanks Even. I don't have access to the machine either as the colleague is moving to another project. I'll have to see if it fails for another *nix user. On 2023-09-24 22:35, Even Rouault wrote: Le 24/09/2023 à 23:22, Jonathan Moules via gdal-dev a écrit : Also check if the environment isn't messed up regarding PROJ and the PROJ_LIB/PROJ_DATA environment variable Thanks Even; sorry, what does this line mean? I'm guessing you're referring to: https://proj.org/en/9.3/usage/environmentvars.html - what would a "messed up" one look like? Hard to say without access to the machine. Perhaps just try to recreate a new Conda env from scratch ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://list/ s.osgeo.org%2Fmailman%2Flistinfo%2Fgdal-dev=05%7C01%7Cjukka.rahko nen%40maanmittauslaitos.fi%7C1f6517888cb34fa40aed08dbc01389dd%7Cc4f8a6 3255804a1c92371d5a571b71fa%7C0%7C0%7C638314962370381350%7CUnknown%7CTW FpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6 Mn0%3D%7C3000%7C%7C%7C=3fVru6K6Ndpkv35FnAbMOlT%2BM96USO7wywqx550 uRUs%3D=0 ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] layer.GetSpatialRef() fails on linux for shapefiles
Well, it seems that PROJ_DATA isn't set in their environment. But it's not set in mine either (`print(os.environ['PROJ_DATA']`)! So no idea why mine works just fine without it (Windows 11 thing?). Creating a PROJ_DATA env var didn't fix anything. Even adding it explicitly in Python. Their log file does have this in at a WARNING level: `PROJ: proj_create_from_database: Open of /home/user/anaconda3/envs/env1/share/proj failed` That path has: `drwxrwxr-x` permissions. To answer Even's earlier question: `ogrinfo /path/to/shape.shp` works fine on their system. On 2023-09-28 12:37, Rahkonen Jukka wrote: Hi, Then they should add that environment if they do not know that they do not belong to "most users" https://proj.org/en/9.3/usage/environmentvars.html -Jukka Rahkonen- -Alkuperäinen viesti- Lähettäjä: gdal-dev Puolesta Jonathan Moules via gdal-dev Lähetetty: torstai 28. syyskuuta 2023 14.10 Vastaanottaja: gdal-dev@lists.osgeo.org; Even Rouault Aihe: Re: [gdal-dev] layer.GetSpatialRef() fails on linux for shapefiles Hi Even, The colleague doesn't have either a PROJ_LIB or a PROJ_DATA environment variable. I asked another colleague to try it; they're on Ubuntu 20.04, and it worked for them. I believe using the same setup instructions. Cheers, Jonathan On 2023-09-24 22:37, Jonathan Moules via gdal-dev wrote: Thanks Even. I don't have access to the machine either as the colleague is moving to another project. I'll have to see if it fails for another *nix user. On 2023-09-24 22:35, Even Rouault wrote: Le 24/09/2023 à 23:22, Jonathan Moules via gdal-dev a écrit : Also check if the environment isn't messed up regarding PROJ and the PROJ_LIB/PROJ_DATA environment variable Thanks Even; sorry, what does this line mean? I'm guessing you're referring to: https://proj.org/en/9.3/usage/environmentvars.html - what would a "messed up" one look like? Hard to say without access to the machine. Perhaps just try to recreate a new Conda env from scratch ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://list/ s.osgeo.org%2Fmailman%2Flistinfo%2Fgdal-dev=05%7C01%7Cjukka.rahko nen%40maanmittauslaitos.fi%7C1f6517888cb34fa40aed08dbc01389dd%7Cc4f8a6 3255804a1c92371d5a571b71fa%7C0%7C0%7C638314962370381350%7CUnknown%7CTW FpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6 Mn0%3D%7C3000%7C%7C%7C=3fVru6K6Ndpkv35FnAbMOlT%2BM96USO7wywqx550 uRUs%3D=0 ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] layer.GetSpatialRef() fails on linux for shapefiles
Hi Even, The colleague doesn't have either a PROJ_LIB or a PROJ_DATA environment variable. I asked another colleague to try it; they're on Ubuntu 20.04, and it worked for them. I believe using the same setup instructions. Cheers, Jonathan On 2023-09-24 22:37, Jonathan Moules via gdal-dev wrote: Thanks Even. I don't have access to the machine either as the colleague is moving to another project. I'll have to see if it fails for another *nix user. On 2023-09-24 22:35, Even Rouault wrote: Le 24/09/2023 à 23:22, Jonathan Moules via gdal-dev a écrit : > Also check if the environment isn't messed up regarding PROJ and the PROJ_LIB/PROJ_DATA environment variable Thanks Even; sorry, what does this line mean? I'm guessing you're referring to: https://proj.org/en/9.3/usage/environmentvars.html - what would a "messed up" one look like? Hard to say without access to the machine. Perhaps just try to recreate a new Conda env from scratch ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] layer.GetSpatialRef() fails on linux for shapefiles
This applies to *all* shapefiles for that user. There are about 10 in the sample and they're from a diverse range of UK organisations. Mostly 27700 (British National Grid), but likely one or two WGS84. Hence thinking it may be environmental, but I don't have any real *nix knowledge to use to diagnose the issue. On 2023-09-24 22:52, Jan Heckman wrote: Sorry to break in, but surely, we would like to see the .prj file in question for a simple try at reproducing? On Sun, Sep 24, 2023 at 11:37 PM Jonathan Moules via gdal-dev wrote: Thanks Even. I don't have access to the machine either as the colleague is moving to another project. I'll have to see if it fails for another *nix user. On 2023-09-24 22:35, Even Rouault wrote: > > Le 24/09/2023 à 23:22, Jonathan Moules via gdal-dev a écrit : >> >> > Also check if the environment isn't messed up regarding PROJ and >> the PROJ_LIB/PROJ_DATA environment variable >> >> Thanks Even; sorry, what does this line mean? I'm guessing you're >> referring to: https://proj.org/en/9.3/usage/environmentvars.html - >> what would a "messed up" one look like? >> > Hard to say without access to the machine. Perhaps just try to > recreate a new Conda env from scratch > > ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] layer.GetSpatialRef() fails on linux for shapefiles
Thanks Even. I don't have access to the machine either as the colleague is moving to another project. I'll have to see if it fails for another *nix user. On 2023-09-24 22:35, Even Rouault wrote: Le 24/09/2023 à 23:22, Jonathan Moules via gdal-dev a écrit : > Also check if the environment isn't messed up regarding PROJ and the PROJ_LIB/PROJ_DATA environment variable Thanks Even; sorry, what does this line mean? I'm guessing you're referring to: https://proj.org/en/9.3/usage/environmentvars.html - what would a "messed up" one look like? Hard to say without access to the machine. Perhaps just try to recreate a new Conda env from scratch ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] layer.GetSpatialRef() fails on linux for shapefiles
> Also check if the environment isn't messed up regarding PROJ and the PROJ_LIB/PROJ_DATA environment variable Thanks Even; sorry, what does this line mean? I'm guessing you're referring to: https://proj.org/en/9.3/usage/environmentvars.html - what would a "messed up" one look like? Thanks, Jonathan On 2023-09-21 13:39, Even Rouault via gdal-dev wrote: Run ogrinfo on one shapefile to see if some more interesting error message is displayed Also check if the environment isn't messed up regarding PROJ and the PROJ_LIB/PROJ_DATA environment variable Le 21/09/2023 à 14:34, Jonathan Moules via gdal-dev a écrit : Yeah, I'm afraid the error message is pretty much non-existent: Traceback (most recent call last): File "/home/name/Code/DSHub/PoC/Extract-Transform-Load-POC/Metadata/AME/src/info_vector.py", line 119, in get_layer_metadata tmp_projection = layer.GetSpatialRef() ^ File "/home/name/anaconda3/envs/AME_env/lib/python3.11/site-packages/osgeo/ogr.py", line 1990, in GetSpatialRef return _ogr.Layer_GetSpatialRef(self, *args) ^ RuntimeError Suggestions welcome. On 18/09/2023 12:51, Javier Jimenez Shaw wrote: Hi Jonathan Which exact RuntimeError are you getting? It can be for several reasons (probably an installation or configuration issue). Best, Javier On Mon, 18 Sept 2023 at 11:06, Jonathan Moules wrote: Hi List, I'm trying to get vector layer information via OGR and Python: ``` layer.GetSpatialRef() ``` This works fine for me on Windows with GDAL 3.7.1 on various different types of files (Shapefile, GPKG, GML, KML, GDB). But for my colleague on Ubuntu 22.0.4.3, also on GDAL 3.7.1 (via Conda), they get a Python RuntimeError for all shapefiles (the exact same shapefiles that work fine for me). It works for Geopackages for them. Anyone have any thoughts? Thanks, Jonathan ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev -- http://www.spatialys.com My software is free, but my time generally not. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] layer.GetSpatialRef() fails on linux for shapefiles
Tested against the following data types (OGR and GDAL): GeoJSON, KML, KMZ, Shapefile, ECW, IMG, TIF, GML, ASC, Geopackage, NetCDF. Only the shapefile fails. On 21/09/2023 13:34, Jonathan Moules wrote: Yeah, I'm afraid the error message is pretty much non-existent: Traceback (most recent call last): File "/home/name/Code/DSHub/PoC/Extract-Transform-Load-POC/Metadata/AME/src/info_vector.py", line 119, in get_layer_metadata tmp_projection = layer.GetSpatialRef() ^ File "/home/name/anaconda3/envs/AME_env/lib/python3.11/site-packages/osgeo/ogr.py", line 1990, in GetSpatialRef return _ogr.Layer_GetSpatialRef(self, *args) ^ RuntimeError Suggestions welcome. On 18/09/2023 12:51, Javier Jimenez Shaw wrote: Hi Jonathan Which exact RuntimeError are you getting? It can be for several reasons (probably an installation or configuration issue). Best, Javier On Mon, 18 Sept 2023 at 11:06, Jonathan Moules wrote: Hi List, I'm trying to get vector layer information via OGR and Python: ``` layer.GetSpatialRef() ``` This works fine for me on Windows with GDAL 3.7.1 on various different types of files (Shapefile, GPKG, GML, KML, GDB). But for my colleague on Ubuntu 22.0.4.3, also on GDAL 3.7.1 (via Conda), they get a Python RuntimeError for all shapefiles (the exact same shapefiles that work fine for me). It works for Geopackages for them. Anyone have any thoughts? Thanks, Jonathan ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] layer.GetSpatialRef() fails on linux for shapefiles
Yeah, I'm afraid the error message is pretty much non-existent: Traceback (most recent call last): File "/home/name/Code/DSHub/PoC/Extract-Transform-Load-POC/Metadata/AME/src/info_vector.py", line 119, in get_layer_metadata tmp_projection = layer.GetSpatialRef() ^ File "/home/name/anaconda3/envs/AME_env/lib/python3.11/site-packages/osgeo/ogr.py", line 1990, in GetSpatialRef return _ogr.Layer_GetSpatialRef(self, *args) ^ RuntimeError Suggestions welcome. On 18/09/2023 12:51, Javier Jimenez Shaw wrote: Hi Jonathan Which exact RuntimeError are you getting? It can be for several reasons (probably an installation or configuration issue). Best, Javier On Mon, 18 Sept 2023 at 11:06, Jonathan Moules wrote: Hi List, I'm trying to get vector layer information via OGR and Python: ``` layer.GetSpatialRef() ``` This works fine for me on Windows with GDAL 3.7.1 on various different types of files (Shapefile, GPKG, GML, KML, GDB). But for my colleague on Ubuntu 22.0.4.3, also on GDAL 3.7.1 (via Conda), they get a Python RuntimeError for all shapefiles (the exact same shapefiles that work fine for me). It works for Geopackages for them. Anyone have any thoughts? Thanks, Jonathan ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] layer.GetSpatialRef() fails on linux for shapefiles
Hi Javier, I don't recall there being any exception message exposed. It was simply a RuntimeError. I guessed it's probably installation/configuration, but surely GDAL should Just Work? Especially if it's coming via Conda? (I don't live in the Linux world; I just used a wheel). Cheers, Jonathan On 18/09/2023 12:51, Javier Jimenez Shaw wrote: Hi Jonathan Which exact RuntimeError are you getting? It can be for several reasons (probably an installation or configuration issue). Best, Javier On Mon, 18 Sept 2023 at 11:06, Jonathan Moules wrote: Hi List, I'm trying to get vector layer information via OGR and Python: ``` layer.GetSpatialRef() ``` This works fine for me on Windows with GDAL 3.7.1 on various different types of files (Shapefile, GPKG, GML, KML, GDB). But for my colleague on Ubuntu 22.0.4.3, also on GDAL 3.7.1 (via Conda), they get a Python RuntimeError for all shapefiles (the exact same shapefiles that work fine for me). It works for Geopackages for them. Anyone have any thoughts? Thanks, Jonathan ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] layer.GetSpatialRef() fails on linux for shapefiles
Hi List, I'm trying to get vector layer information via OGR and Python: ``` layer.GetSpatialRef() ``` This works fine for me on Windows with GDAL 3.7.1 on various different types of files (Shapefile, GPKG, GML, KML, GDB). But for my colleague on Ubuntu 22.0.4.3, also on GDAL 3.7.1 (via Conda), they get a Python RuntimeError for all shapefiles (the exact same shapefiles that work fine for me). It works for Geopackages for them. Anyone have any thoughts? Thanks, Jonathan ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Corrupted geopackage
Hi Isabel,
If you don't get an answer here, you may want to try the SQLite forum:
https://sqlite.org/forum/forum - If anyone knows, it'll be someone there.
Maybe also of interest: https://sqlite.org/howtocorrupt.html#cfgerr
Cheers,
Jonathan
On 2022-02-14 08:17, Isabel Kiefer wrote:
Hi everyone,
I've a problem with a corrupted geopackage.
It should contain a table with 694 entries, but only 78 are visible
when opening the file with DB Browser for SQLite or QGIS. The table
gpkg_ogr_contents says 694 though. When opening the .gpkg with a
Notepad or similar, I can see that there are more than 78 entries.
The error message in DB Browser is "database disk image is malformed
in "PRAGMA "main".TABLE_INFO("rtree_XXX");" So it seems that there is
a problem with the index of the concerned table.
Does anyone know how to fix this? I tried
https://stackoverflow.com/questions/18259692/how-to-recover-a-corrupt-sqlite3-database/18260642
but
without success.
Thanks in advance for your help!
Isabel
--
Isabel Kiefer
OPENGIS.ch
isa...@opengis.ch
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/gdal-dev___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Reading a JP2 file - obj_decode() failed
I guess this is what you want: 132000 248000 Which I'm assuming means it's not tile. On 2020-08-09 15:46, Even Rouault wrote: On dimanche 9 août 2020 15:30:12 CEST Jonathan Moules wrote: > Hi Even, > > Ah, it's comma delimited. I tried spaces. > > Here's the response, with with light editing to get rid of duplicate > info (The bands are all the same). Nothing in the output with the string > "tile". Ah, sorry, I now see those debug traces aren't compiled in in release mode. OK, so download https://github.com/OSGeo/gdal/blob/master/gdal/swig/python/samples/dump_jp2.py and run python dump_jp2.py FILENAME.jp2 | grep Tsiz -- Spatialys - Geospatial professional services http://www.spatialys.com ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Reading a JP2 file - obj_decode() failed
Hi Even, Ah, it's comma delimited. I tried spaces. Here's the response, with with light editing to get rid of duplicate info (The bands are all the same). Nothing in the output with the string "tile". Cheers, Jonathan OPENJPEG: info: Start to read j2k main header (0). OPENJPEG: info: Main header has been correctly decoded. OPENJPEG: SRGB color space GDALJP2Metadata: Got projection from GeoJP2 (geotiff) box (0): PROJCS["British National Grid (ORD SURV GB)",GEOGCS["OSGB 1936",DATUM["OSGB_1936",SPHEROID["Airy1830",6377563.396,299.3249612664951,AUTHORITY["EPSG","7001"]],AUTHORITY["EPSG","6277"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4277"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",49],PARAMETER["central_meridian",-2],PARAMETER["scale_factor",0.9996012717],PARAMETER["false_easting",40],PARAMETER["false_northing",-10],UNIT["metre",1,AUTHORITY["EPSG","9001"]]] MDReaderPleiades: Not a Pleiades product MDReaderPleiades: Not a Pleiades product GDAL: GDALOpen(FILENAME.jp2, this=060407F0) succeeds as JP2OpenJPEG. Driver: JP2OpenJPEG/JPEG-2000 driver based on OpenJPEG library GDAL: GDALDefaultOverviews::OverviewScan() Files: FILENAME.jp2 FILENAME.jp2.aux.xml ... Metadata: GEOTIFF_CHAR__GeogAngularUnitsGeoKey=Angular_Degree GEOTIFF_CHAR__GeogEllipsoidGeoKey=Ellipse_Airy_1830 GEOTIFF_CHAR__GeogGeodeticDatumGeoKey=Datum_OSGB_1936 GEOTIFF_CHAR__GeographicTypeGeoKey=GCS_OSGB_1936 GEOTIFF_CHAR__GTModelTypeGeoKey=ModelTypeProjected GEOTIFF_CHAR__GTRasterTypeGeoKey=RasterPixelIsArea GEOTIFF_CHAR__ProjCoordTransGeoKey=CT_TransverseMercator GEOTIFF_CHAR__ProjectedCSTypeGeoKey=User-Defined GEOTIFF_CHAR__ProjectionGeoKey=User-Defined GEOTIFF_CHAR__ProjLinearUnitsGeoKey=Linear_Meter GEOTIFF_NUM__1024__GTModelTypeGeoKey=1 GEOTIFF_NUM__1025__GTRasterTypeGeoKey=1 GEOTIFF_NUM__1026__GTCitationGeoKey=British National Grid (ORD SURV GB) GEOTIFF_NUM__2048__GeographicTypeGeoKey=4277 GEOTIFF_NUM__2049__GeogCitationGeoKey=OSGB 1936 GEOTIFF_NUM__2050__GeogGeodeticDatumGeoKey=6277 GEOTIFF_NUM__2054__GeogAngularUnitsGeoKey=9102 GEOTIFF_NUM__2056__GeogEllipsoidGeoKey=7001 GEOTIFF_NUM__2057__GeogSemiMajorAxisGeoKey=6377563.396000 GEOTIFF_NUM__2059__GeogInvFlatteningGeoKey=299.324961 GEOTIFF_NUM__3072__ProjectedCSTypeGeoKey=32767 GEOTIFF_NUM__3074__ProjectionGeoKey=32767 GEOTIFF_NUM__3075__ProjCoordTransGeoKey=1 GEOTIFF_NUM__3076__ProjLinearUnitsGeoKey=9001 GEOTIFF_NUM__3080__ProjNatOriginLongGeoKey=-2.00 GEOTIFF_NUM__3081__ProjNatOriginLatGeoKey=49.00 GEOTIFF_NUM__3082__ProjFalseEastingGeoKey=40.00 GEOTIFF_NUM__3083__ProjFalseNorthingGeoKey=-10.00 GEOTIFF_NUM__3092__ProjScaleAtNatOriginGeoKey=0.999601 Image Structure Metadata: INTERLEAVE=PIXEL OGRCT: Source: +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=40 +y_0= -10 +datum=OSGB36 +units=m +no_defs OGRCT: Target: +proj=longlat +datum=OSGB36 +no_defs Corner Coordinates: ... Band 1 Block=1024x1024 Type=Byte, ColorInterp=Red Min=168.000 Max=255.000 Minimum=168.000, Maximum=255.000, Mean=243.359, StdDev=16.164 Overviews: 66000x124000, 33000x62000, 16500x31000, 8250x15500, 4125x7750, 2062 x3875, 1031x1937, 515x968, 257x484, 128x242, 64x121 Overviews: arbitrary Metadata: STATISTICS_APPROXIMATE=YES STATISTICS_MAXIMUM=255 STATISTICS_MEAN=243.35889029004 STATISTICS_MINIMUM=168 STATISTICS_STDDEV=16.164314681862 STATISTICS_VALID_PERCENT=100 Image Structure Metadata: COMPRESSION=JPEG2000 ... GDAL: GDALClose(FILENAME.jp2,this=060407F0) On 2020-08-09 15:21, Even Rouault wrote: On dimanche 9 août 2020 12:47:33 CEST Jonathan Moules wrote: > Hi Even, > Thanks for getting back to me. I'm a little rusty with GDAL; it has > probably been 5+ years since I last used it. > > I can confirm my QGIS install version (OSGEO4W) does have the JP2ECW > driver (`gdalinfo --formats`). > > I'm using the Python GDAL build from > https://www.lfd.uci.edu/~gohlke/pythonlibs/#gdal - doesn't say it > includes JP2ECW. I only have access to that GDAL build via Python but > haven't managed to figure out what the Python API equivalent of > "gdalinfo --formats" is (assuming there is one). > > I don't seem to be able to force JP2OpenJPEG on the OSGeo4W one. > gdalinfo --config GDAL_SKIP JP2ECW FILENAME.jp2 = uses the MrSID driver > gdalinfo --config GDAL_SKIP JP2ECW --config GDAL_SKIP JP2MrSID > FILENAME.jp2 = Goes back to JP2ECW If you specify several times a configuration option, the last value for a given key will override the previous one(s). You can skip several drivers by using
Re: [gdal-dev] Reading a JP2 file - obj_decode() failed
Hi Even,
Thanks for getting back to me. I'm a little rusty with GDAL; it has
probably been 5+ years since I last used it.
I can confirm my QGIS install version (OSGEO4W) does have the JP2ECW
driver (`gdalinfo --formats`).
I'm using the Python GDAL build from
https://www.lfd.uci.edu/~gohlke/pythonlibs/#gdal - doesn't say it
includes JP2ECW. I only have access to that GDAL build via Python but
haven't managed to figure out what the Python API equivalent of
"gdalinfo --formats" is (assuming there is one).
I don't seem to be able to force JP2OpenJPEG on the OSGeo4W one.
gdalinfo --config GDAL_SKIP JP2ECW FILENAME.jp2 = uses the MrSID driver
gdalinfo --config GDAL_SKIP JP2ECW --config GDAL_SKIP JP2MrSID
FILENAME.jp2 = Goes back to JP2ECW
(Doesn't support the -if flag).
But I found `gdal.Info("FILENAME.jp2")` on the Python version, and that
does use the JP2 OpenJPEG driver. No tile information on there:
Driver: JP2OpenJPEG/JPEG-2000 driver based on OpenJPEG library
...
Image Structure Metadata:
INTERLEAVE=PIXEL
...
Band 1 Block=1024x1024 Type=Byte, ColorInterp=Red
Overviews: 66000x124000, 33000x62000, 16500x31000, 8250x15500,
4125x7750, 2063x3875, 1032x1938, 516x969, 258x485, 129x243, 65x122
Overviews: arbitrary
Image Structure Metadata:
COMPRESSION=JPEG2000
...
I can't help but notice the block size and overview levels are reporting
different info compared to the JP2ECW output (ECW one has different
rounding <= 2062*3875):
Band 1 Block=256x256 Type=Byte, ColorInterp=Red
Description = Red
Overviews: 66000x124000, 33000x62000, 16500x31000, 8250x15500,
4125x7750, 2062
x3875, 1031x1937, 515x968, 257x484
and JP2MrSID for comparison:
Band 1 Block=1024x128 Type=Byte, ColorInterp=Red
Min=168.000 Max=255.000
Minimum=168.000, Maximum=255.000, Mean=243.359, StdDev=16.164
Overviews: 66000x124000, 33000x62000, 16500x31000, 8250x15500,
4125x7750, 2063
x3875, 1032x1938, 516x969, 258x485, 129x243, 65x122, 33x61
Metadata:
STATISTICS_APPROXIMATE=YES
STATISTICS_MAXIMUM=255
STATISTICS_MEAN=243.35889029004
STATISTICS_MINIMUM=168
STATISTICS_STDDEV=16.164314681862
STATISTICS_VALID_PERCENT=100
All three report different Overview levels too. Odd but no idea if it is
pertinent to this problem.
Cheers,
Jonathan
p.s. I was looking for Python API docs but it seems there aren't any.
This page:
https://pypi.org/project/GDAL/#usage - has a link to
http://www.gdal.org/gdal_tutorial.html, but that's a 404. I think it
wants to be: https://gdal.org/tutorials/index.html
On 2020-08-09 11:14, Even Rouault wrote:
Jonathan,
The GDAL build that comes with rasterio must have only JP2OpenJPEG
enabled (the error with opj_decode() comes from that driver), whereas
the GDAL build in QGIS has also JP2ECW, which is used in priority.
The file is large. It would be interesting to know if it is tiled. If
you do "gdalinfo --format FILENAME.jp2 --config GDAL_SKIP JP2ECW", and
look at the beginning of the output, you should see lines like
OPENJPEG: nTileW =
OPENJPEG: nTileH =
which values are reported ?
If it is tiled (values <= 2048 let's say), this should normally work.
Openjpeg will have more difficulties with single-tile JPEG2000 files,
although latest versions (2.3.1) have seen improvements in that regard
Even
--
Spatialys - Geospatial professional services
http://www.spatialys.com
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] Reading a JP2 file - obj_decode() failed
Hi List, I have a JP2 file I'm trying to process with Rasterio in Python. But I'm getting this error: rasterio.errors.RasterioIOError: Read or write failed. FILENAME.jp2, band 1: IReadBlock failed at X offset 40, Y offset 90: opj_decode() failed That error seems to come from GDAL. I've tried setting the driver to both "JPEG2000" and "JP2OpenJPEG" to no avail. The file opens fine in QGIS, and I tried both gdal2tiles and rgb2pct (via the QGIS installation, which is different from the python/rasterio install): those both worked fine too, so clearly GDAL can read it when standalone. Any suggestions? Thanks, Jonathan p.s. - Also tried the JP2ECW driver but couldn't figure out how to link it to GDAL. I mention it only because the links to the JP2ECW on https://gdal.org/drivers/raster/jp2ecw.html are now invalid. They seem to have moved stuff without using redirects. Seems to be here now: https://download.hexagongeospatial.com/search?keyword==en_family=29e5f4b2-f943-4710-9b4a-bd1d16525d1e_version=a0edddba-0ce6-4cfb-bbe4-95e59d3c87d2 gdalinfo output below: Driver: JP2ECW/ERDAS JPEG2000 (SDK 5.3) Files: FILENAME.jp2 Size is 132000, 248000 Coordinate System is: PROJCS["British National Grid (ORD SURV GB)", GEOGCS["OSGB 1936", DATUM["OSGB_1936", SPHEROID["Airy 1830",6377563.396,299.3249612664951, AUTHORITY["EPSG","7001"]], AUTHORITY["EPSG","6277"]], PRIMEM["Greenwich",0], UNIT["degree",0.0174532925199433], AUTHORITY["EPSG","4277"]], PROJECTION["Transverse_Mercator"], PARAMETER["latitude_of_origin",49], PARAMETER["central_meridian",-2], PARAMETER["scale_factor",0.9996012717], PARAMETER["false_easting",40], PARAMETER["false_northing",-10], UNIT["metre",1, AUTHORITY["EPSG","9001"]]] Origin = (0.000,124.000) Pixel Size = (5.000,-5.000) Metadata: COLORSPACE=RGB COMPRESSION_RATE_TARGET=0 Corner Coordinates: Upper Left ( 0.000, 124.000) ( 9d22'13.04"W, 60d50'27.04"N) Lower Left ( 0.000, 0.000) ( 7d33'23.21"W, 49d45'58.27"N) Upper Right ( 66.000, 124.000) ( 2d48'15.07"E, 60d57'26.43"N) Lower Right ( 66.000, 0.000) ( 1d37' 0.77"E, 49d50'35.22"N) Center ( 33.000, 62.000) ( 3d 6'26.52"W, 55d28' 7.68"N) Band 1 Block=256x256 Type=Byte, ColorInterp=Red Description = Red Overviews: 66000x124000, 33000x62000, 16500x31000, 8250x15500, 4125x7750, 2062 x3875, 1031x1937, 515x968, 257x484 Band 2 Block=256x256 Type=Byte, ColorInterp=Green Description = Green Overviews: 66000x124000, 33000x62000, 16500x31000, 8250x15500, 4125x7750, 2062 x3875, 1031x1937, 515x968, 257x484 Band 3 Block=256x256 Type=Byte, ColorInterp=Blue Description = Blue Overviews: 66000x124000, 33000x62000, 16500x31000, 8250x15500, 4125x7750, 2062 x3875, 1031x1937, 515x968, 257x484 ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Motion: adopt RFC 72: Run tests with pytest
Hi, PyTest is a great test-suite. If I may make one suggestion as someone who has used it for a while - it's worth spending a little thought on coming up with a scheme for test-ids. Especially if you're going to use parameterisation. Otherwise PyTest comes up with names that may be accurate (they're a concatenation of the parameters), but are relatively meaningless. For example "gdaladdo--100", "gdaladdo-foo", etc, as compared to more useful ids like "gdaladdo-too-low-input", "gdaladdo-bad-string-input" which tell you immediately what the test is actually meant to be testing. There's a hook called pytest_make_parametrize_id which allows you to create your own ids (I find the built-in methods of id generation either cumbersome, or bit-rot prone). I'm suggesting it now because it's more helpful if you do it from the start. Cheers, Jonathan On 2018-12-05 22:40, Craig de Stigter wrote: Hi I appreciate your comments on the pytest proposal and all the support to help get it this far. Given no actionable improvements have been suggested, and the feedback thus far seems encouraging... I move to adopt RFC 72: Run tests with pytest. https://trac.osgeo.org/gdal/wiki/rfc72_pytest Cheers Craig de Stigter ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] ElasticSearch and native geometry
Hi Jukka, I'm still very new to ElasticSearch (so take with plenty of salt!), but if you want to keep the geom in _source and calculate on the fly, I guess you could maybe try it with one of the scripting languages: https://www.elastic.co/guide/en/elasticsearch/reference/current/modules-scripting.html - Python's on there so then you'd have access to all of the py-geo stuff. Of course, assuming it worked, the trade-off would be reduced storage for (considerably?) increase processing usage and complexity. The other option would be to store the native geometry as a "binary" datatype - https://www.elastic.co/guide/en/elasticsearch/reference/current/binary.html Cheers, Jonathan On 23/11/2016 15:41, Rahkonen Jukka (MML) wrote: Hi, What if somebody would like to use ElasticSearch for queries but still somehow get the native geometries included in the output? The geometry must be re-projected into EPSG:4326 for creating either geo_point or geo_shape, but can anybody suggest a clever way for saving the native geometry? It is of course possible to save it into a text field for example as WKT of JSON but could it be possible to keep the original geometry in “_source” and just create the geo_point or geo_shape on-the-fly or the index? -Jukka Rahkonen- ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] Verbose mode in GDAL
Hi List, I'm running the following four queries against the same data (a packbits compressed geotiff): gdaladdo -r average --config COMPRESS_OVERVIEW JPEG --config PHOTOMETRIC_OVERVIEW YCBCR i1.tif 2 4 8 16 gdaladdo -r average --config PHOTOMETRIC_OVERVIEW YCBCR i2.tif 2 4 8 16 gdaladdo -r average i3.tif 2 4 8 16 gdaladdo -r average --config COMPRESS_OVERVIEW PACKBITS i4.tif 2 4 8 16 All four outputs are binary identical and contain the requested internal overviews. I can think of several possible reasons for this outcome (different compression type seems most probable), but invisibly failing to action semantically valid flags isn't good for raising or debugging errors. Does GDAL have a "verbose" mode? I can't seem to find one. Cheers, Jonathan ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Standard Parallel 1 - being changed but not displayed?
Hi Even,Ok, that clears things up a little. I guess the real question I have at
this point is - which is correct in this circumstance?
The original GeoTIFF?
ArcGIS?
GDAL's conversion?
None of the above? I get the feeling it's all a bit subjective.
>From your comments, I'm surprised this doesn't happen more often with various
>packages\files. I thought projections were nice and reasonably well defined
>things.
Thanks again,
Jonathan
On Tue, 17 Nov 2015 10:33:37 -0800 Even
Rouaulteven.roua...@spatialys.com wrote
Le mardi 17 novembre 2015 19:13:52, Jonathan Moules a écrit :
Hi Even,
Thanks for the reply.
I guess GDAL chose Mercator_1SP because that's the one that's explicitly
defined in the projection text.
No, it is not explicitly defined, and GDAL ignores PCSCitationGeoKey to
determine the projection method. It will use ProjCoordTransGeoKey =
CT_Mercator for that, and in GDAL 1.11, it can only map it to Mercator_1SP
I've had a quick look at the two files using listgeo, and highlighted the
different lines with prefix below (it seems GDAL removes 6 tags, and
adds one).
The before:
Keyed_Information:
GTModelTypeGeoKey (Short,1): ModelTypeProjected
GTRasterTypeGeoKey (Short,1): RasterPixelIsArea
GTCitationGeoKey (Ascii,9): "Mercator"
GeographicTypeGeoKey (Short,1): GCS_WGS_84
GeogCitationGeoKey (Ascii,7): "WGS 84"
GeogGeodeticDatumGeoKey (Short,1): Datum_WGS84
GeogLinearUnitsGeoKey (Short,1): Linear_Meter
GeogAngularUnitsGeoKey (Short,1): Angular_Degree
GeogEllipsoidGeoKey (Short,1): Ellipse_WGS_84
GeogSemiMajorAxisGeoKey (Double,1): 6378137
GeogSemiMinorAxisGeoKey (Double,1): 6356752.31424518
GeogInvFlatteningGeoKey (Double,1): 298.257223563
ProjectedCSTypeGeoKey (Short,1): User-Defined
PCSCitationGeoKey (Ascii,9): "Mercator"
ProjectionGeoKey (Short,1): User-Defined
ProjCoordTransGeoKey (Short,1): CT_Mercator
ProjLinearUnitsGeoKey (Short,1): Linear_Meter
ProjStdParallel1GeoKey (Double,1): 60
ProjNatOriginLongGeoKey (Double,1): 0
ProjNatOriginLatGeoKey (Double,1): 0
ProjFalseEastingGeoKey (Double,1): 0
ProjFalseNorthingGeoKey (Double,1): 0
End_Of_Keys.
The after:
Keyed_Information:
GTModelTypeGeoKey (Short,1): ModelTypeProjected
GTRasterTypeGeoKey (Short,1): RasterPixelIsArea
GTCitationGeoKey (Ascii,9): "Mercator"
GeographicTypeGeoKey (Short,1): GCS_WGS_84
GeogCitationGeoKey (Ascii,7): "WGS 84"
GeogAngularUnitsGeoKey (Short,1): Angular_Degree
GeogSemiMajorAxisGeoKey (Double,1): 6378137
GeogInvFlatteningGeoKey (Double,1): 298.257223563
ProjectedCSTypeGeoKey (Short,1): User-Defined
ProjectionGeoKey (Short,1): User-Defined
ProjCoordTransGeoKey (Short,1): CT_Mercator
ProjLinearUnitsGeoKey (Short,1): Linear_Meter
ProjNatOriginLongGeoKey (Double,1): 0
ProjNatOriginLatGeoKey (Double,1): 0
ProjFalseEastingGeoKey (Double,1): 0
ProjFalseNorthingGeoKey (Double,1): 0
ProjScaleAtNatOriginGeoKey (Double,1): 1
End_Of_Keys.
End_Of_Geotiff.
Actually when trying that with GDAL 2.0, you would now get :
PROJCS["Mercator",
GEOGCS["WGS 84",
DATUM["WGS_1984",
SPHEROID["WGS 84",6378137,298.257223563,
AUTHORITY["EPSG","7030"]],
AUTHORITY["EPSG","6326"]],
PRIMEM["Greenwich",0],
UNIT["degree",0.0174532925199433],
AUTHORITY["EPSG","4326"]],
PROJECTION["Mercator_2SP"],
PARAMETER["standard_parallel_1",60],
PARAMETER["central_meridian",0],
PARAMETER["false_easting",0],
PARAMETER["false_northing",0],
UNIT["metre",1,
AUTHORITY["EPSG","9001"]]]
due to the changes done in https://trac.osgeo.org/gdal/ticket/5791
The new logic to determine Mercator_1SP vs _2SP is best explainted by the code
itself ;-)
{{{
/* If a lat_ts was specified use 2SP, otherwise use 1SP */
if (psDefn-ProjParmId[2] == ProjStdParallel1GeoKey)
{
if (psDefn-ProjParmId[4] == ProjScaleAtNatOriginGeoKey)
CPLError( CE_Warning, CPLE_AppDefined,
"Mercator projection should not define both
StdParallel1 and ScaleAtNatOrigin.\n"
"Using StdParallel1 and ignoring
ScaleAtNatOrigin.\n" );
oSRS.SetMercator2SP( adfParm[2],
adfParm[0], adfParm[1],
adfParm[5], adfParm[6]);
}
else
oSRS.SetMercator( adfParm[0], adfParm[1],
adfParm[4],
adfParm[5], adfParm[6] );
}}}
That new behaviour is more consistant with the geotiff keys (and should be the
one you get with ArcGIS initialy), although it perhaps reveals a problem with
the encoding of the geotiff file itself if you say that the Mercator_1SP
interpretation with scale=1 leads to correct geopositioning.
I'm not clear on why the things that changed were changed, but I can see
that
[gdal-dev] Standard Parallel 1 - being changed but not displayed?
Hi List, I have a Geotiff which includes this projection: PROJ.4 : '+proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs ' OGC WKT : PROJCS["Mercator", GEOGCS["WGS 84", DATUM["WGS_1984", SPHEROID["WGS 84",6378137,298.257223563, AUTHORITY["EPSG","7030"]], AUTHORITY["EPSG","6326"]], PRIMEM["Greenwich",0], UNIT["degree",0.0174532925199433], AUTHORITY["EPSG","4326"]], PROJECTION["Mercator_1SP"], PARAMETER["central_meridian",0], PARAMETER["scale_factor",1], PARAMETER["false_easting",0], PARAMETER["false_northing",0], UNIT["metre",1, AUTHORITY["EPSG","9001"]]] If I load this raster into ArcGIS, it displays in the wrong place (a few thousand kilometres North). I then run it through gdal_translate (GDAL 1.11.1), with no flags: * gdal_trainslate input.tif output.tif For output.tif, GDALSRSInfo shows that the projection is identical, but now the file loads correctly in ArcGIS. The same file works fine in QGIS both before and after the "translation". Looking at the projection info in ArcGIS, it displays one difference: Before (not working): Standard_parallel_1 = 60 After (working): Standard_parallel_1 = 0 But I don't see anything about those in either of the GDALSRSInfo outputs. So my questions: - What is gdal_translate doing to the file to "fix" it? - If it is something to do with Standard Parallel 1 - why isn't this component of the projection exposed by GDAL? Thoughts welcome. Thanks, Jonathan ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Standard Parallel 1 - being changed but not displayed?
Hi Even, Thanks for the reply. I guess GDAL chose Mercator_1SP because that's the one that's explicitly defined in the projection text. I've had a quick look at the two files using listgeo, and highlighted the different lines with prefix below (it seems GDAL removes 6 tags, and adds one). The before: Keyed_Information: GTModelTypeGeoKey (Short,1): ModelTypeProjected GTRasterTypeGeoKey (Short,1): RasterPixelIsArea GTCitationGeoKey (Ascii,9): "Mercator" GeographicTypeGeoKey (Short,1): GCS_WGS_84 GeogCitationGeoKey (Ascii,7): "WGS 84" GeogGeodeticDatumGeoKey (Short,1): Datum_WGS84 GeogLinearUnitsGeoKey (Short,1): Linear_Meter GeogAngularUnitsGeoKey (Short,1): Angular_Degree GeogEllipsoidGeoKey (Short,1): Ellipse_WGS_84 GeogSemiMajorAxisGeoKey (Double,1): 6378137 GeogSemiMinorAxisGeoKey (Double,1): 6356752.31424518 GeogInvFlatteningGeoKey (Double,1): 298.257223563 ProjectedCSTypeGeoKey (Short,1): User-Defined PCSCitationGeoKey (Ascii,9): "Mercator" ProjectionGeoKey (Short,1): User-Defined ProjCoordTransGeoKey (Short,1): CT_Mercator ProjLinearUnitsGeoKey (Short,1): Linear_Meter ProjStdParallel1GeoKey (Double,1): 60 ProjNatOriginLongGeoKey (Double,1): 0 ProjNatOriginLatGeoKey (Double,1): 0 ProjFalseEastingGeoKey (Double,1): 0 ProjFalseNorthingGeoKey (Double,1): 0 End_Of_Keys. The after: Keyed_Information: GTModelTypeGeoKey (Short,1): ModelTypeProjected GTRasterTypeGeoKey (Short,1): RasterPixelIsArea GTCitationGeoKey (Ascii,9): "Mercator" GeographicTypeGeoKey (Short,1): GCS_WGS_84 GeogCitationGeoKey (Ascii,7): "WGS 84" GeogAngularUnitsGeoKey (Short,1): Angular_Degree GeogSemiMajorAxisGeoKey (Double,1): 6378137 GeogInvFlatteningGeoKey (Double,1): 298.257223563 ProjectedCSTypeGeoKey (Short,1): User-Defined ProjectionGeoKey (Short,1): User-Defined ProjCoordTransGeoKey (Short,1): CT_Mercator ProjLinearUnitsGeoKey (Short,1): Linear_Meter ProjNatOriginLongGeoKey (Double,1): 0 ProjNatOriginLatGeoKey (Double,1): 0 ProjFalseEastingGeoKey (Double,1): 0 ProjFalseNorthingGeoKey (Double,1): 0 ProjScaleAtNatOriginGeoKey (Double,1): 1 End_Of_Keys. End_Of_Geotiff. I'm not clear on why the things that changed were changed, but I can see that GDAL removed the ProjStdParallel1GeoKey and value. Should not that have been displayed by GDALSRSInfo for the Before file, even if it was wrongly set for this given projection? Thanks, Jonathan On Tue, 17 Nov 2015 09:18:49 -0800 Even Rouaulteven.roua...@spatialys.com wrote ---- Le mardi 17 novembre 2015 17:55:15, Jonathan Moules a écrit : Hi List, I have a Geotiff which includes this projection: PROJ.4 : '+proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs ' OGC WKT : PROJCS["Mercator", GEOGCS["WGS 84", DATUM["WGS_1984", SPHEROID["WGS 84",6378137,298.257223563, AUTHORITY["EPSG","7030"]], AUTHORITY["EPSG","6326"]], PRIMEM["Greenwich",0], UNIT["degree",0.0174532925199433], AUTHORITY["EPSG","4326"]], PROJECTION["Mercator_1SP"], PARAMETER["central_meridian",0], PARAMETER["scale_factor",1], PARAMETER["false_easting",0], PARAMETER["false_northing",0], UNIT["metre",1, AUTHORITY["EPSG","9001"]]] If I load this raster into ArcGIS, it displays in the wrong place (a few thousand kilometres North). I then run it through gdal_translate (GDAL 1.11.1), with no flags: * gdal_trainslate input.tif output.tif For output.tif, GDALSRSInfo shows that the projection is identical, but now the file loads correctly in ArcGIS. The same file works fine in QGIS both before and after the "translation". Looking at the projection info in ArcGIS, it displays one difference: Before (not working): Standard_parallel_1 = 60 After (working): Standard_parallel_1 = 0 But I don't see anything about those in either of the GDALSRSInfo outputs. So my questions: - What is gdal_translate doing to the file to "fix" it? - If it is something to do with Standard Parallel 1 - why isn't this component of the projection exposed by GDAL? Yes, in Mercator_1SP, there's no Standard Parallel 1, this is for Mercator_2SP. See http://www.remotesensing.org/geotiff/proj_list/mercator_1sp.html http://www.remotesensing.org/geotiff/proj_list/mercator_2sp.html I guess your original geotiff file has some unusual formulation which is detected as Mercator_1SP by GDAL, and probably Mercator_2SP by ArCGIS. You could try with the listgeo ut
Re: [gdal-dev] [Qgis-user] Generating GeoTIFF overview (gdaladdo -r gauss ...) takes days for 2GB file...
Hi Stefan, How large is the image? As in, pixel width height. Also the GDALinfo would be helpful. 7 days is a very long time. I've seen times like that myself for images that were 100,000*50,000 or similar dimensions where there was (still is?) an issue with the pyramiding function in GDAL (it was re-reading the data again and again). Cheers, Jonathan -Original Message- From: qgis-user-boun...@lists.osgeo.org [mailto:qgis-user-boun...@lists.osgeo.org] On Behalf Of Stefan Keller Sent: Sunday, October 12, 2014 1:34 PM To: gdal-dev@lists.osgeo.org; qgis-user Subject: [Qgis-user] Generating GeoTIFF overview (gdaladdo -r gauss ...) takes days for 2GB file... Hi, I have a GeoTIFF file of 2GB size. Now, this question is about gdaladdo and how to optimize loading time of such a large raster file. According to my observations QGIS has fastest loading times when it's a single GeoTIFF file with overview images. Q.1: Does anybody have tipps for even faster solutions loading raster files into in QGIS (prepared by OGR or other preprocessors)? The command gdaladdo -r gauss C:\SG.tif 2 4 8 16 32 64 128 256 takes actually about 7 days on a Lenovo notebook. Q2.: Any hint's on how to accelerate this (or to other tools generating overviews/pyramids)? Yours, S. ___ Qgis-user mailing list qgis-u...@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/qgis-user This message has been scanned for viruses by MailControl - www.mailcontrol.com Click https://www.mailcontrol.com/sr/pdM3aunlCjzGX2PQPOmvUgEBY15Clgt1cLVtYJoNw+TFapmLjYzc0LQlVsaeWay487vQIAp6KUWSPFO0UaI8cQ== to report this email as spam. HR Wallingford and its subsidiaries uses faxes and emails for confidential and legally privileged business communications. They do not of themselves create legal commitments. Disclosure to parties other than addressees requires our specific consent. We are not liable for unauthorised disclosures nor reliance upon them. If you have received this message in error please advise us immediately and destroy all copies of it. HR Wallingford Limited Howbery Park, Wallingford, Oxfordshire, OX10 8BA, United Kingdom Registered in England No. 02562099 ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] Questions about processing four band (RGBI) imagery
Hi List, I have four band (RGB plus Infrared) imagery from some aerial photography we've had flown. I've already used GDAL to create suitable optimised RGB GeoTIFFs using this data and now want to do it with the Infrared band too. The RGBI data is about 1TB uncompressed, currently stored as 2732 GeoTIFF tiles compressed with LZW: Files: SK1400.tif SK1400.tfw Size is 8000, 8000 Coordinate System is: PROJCS[OSGB 1936 / British National Grid, GEOGCS[OSGB 1936, DATUM[OSGB_1936, SPHEROID[Airy 1830,6377563.396,299.3249753150316, AUTHORITY[EPSG,7001]], AUTHORITY[EPSG,6277]], PRIMEM[Greenwich,0], UNIT[degree,0.0174532925199433], AUTHORITY[EPSG,4277]], PROJECTION[Transverse_Mercator], PARAMETER[latitude_of_origin,49], PARAMETER[central_meridian,-2], PARAMETER[scale_factor,0.9996012717], PARAMETER[false_easting,40], PARAMETER[false_northing,-10], UNIT[metre,1, AUTHORITY[EPSG,9001]], AUTHORITY[EPSG,27700]] Origin = (414000.000,301000.000) Pixel Size = (0.125,-0.125) Metadata: AREA_OR_POINT=Area Image Structure Metadata: COMPRESSION=LZW INTERLEAVE=PIXEL Corner Coordinates: Upper Left ( 414000.000, 301000.000) ( 1d47'35.68W, 52d36'22.49N) Lower Left ( 414000.000, 30.000) ( 1d47'35.83W, 52d35'50.12N) Upper Right ( 415000.000, 301000.000) ( 1d46'42.51W, 52d36'22.40N) Lower Right ( 415000.000, 30.000) ( 1d46'42.68W, 52d35'50.03N) Center ( 414500.000, 300500.000) ( 1d47' 9.18W, 52d36' 6.26N) Band 1 Block=8000x32 Type=Byte, ColorInterp=Red Band 2 Block=8000x32 Type=Byte, ColorInterp=Green Band 3 Block=8000x32 Type=Byte, ColorInterp=Blue Band 4 Block=8000x32 Type=Byte, ColorInterp=Undefined My goal is to create a mosaiced GeoTIFF image with all four bands in it. It should be as small as possible without serious quality degredation. I've got some some questions as a result of my experiments so far: - -co PHOTOMETRIC=YCBCR - in gdal_translate this will fail saying it only works with 3 bands. Yet with gdaladdo, compressing overviews with --config PHOTOMETRIC_OVERVIEW YCBCR works fine (or at least gives no errors). Any reason for this? - Is there any special compression to use for four band imagery as I can't use YCBCR? - What's the best way to process the data? I'm using: gdal_translate aerial.vrt aerial.tif -of GTiff -co BIGTIFF=YES -co TILED=YES -co blockxsize=512 -co blockysize=512 -co COMPRESS=JPEG -co JPEG_QUALITY=50 -a_srs EPSG:27700 gdaladdo aerial.tif -r average --config COMPRESS_OVERVIEW JPEG --config JPEG_QUALITY_OVERVIEW 50 --config INTERLEAVE_OVERVIEW PIXEL --config PHOTOMETRIC_OVERVIEW YCBCR 2 4 8 16 32 - How do I get GDAL to *not *set Band 4 to Alpha (which is what the above does). Using GDAL 1.10.1 Thanks, Jonathan -- This transmission is intended for the named addressee(s) only and may contain confidential, sensitive or personal information and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Gdal GeoTIFFs in MapInfo
Hi Even, Well, I don't know MapInfo limits with raster files. Perhaps is it just big raster dimensions that cause problem, or tiling, or .. ? I guess it's size. I've tried with no compression, no tiles, and both no compression or tiling, all fail. So I guess a MapInfo thing then. Thanks, Jonathan -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Extra Decimal places added to ASCII gdal_translate mosaic
Hi David, Here's a crazy idea could you perhaps use mm and then output to an integer raster? No problem representing 50597mm. I like you thinking! But on reflection suspect that would confound a few too many of our users. It's certainly a nice lateral-thinking solution to bear in mind for the future though. Thanks! Jonathan -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Extra Decimal places added to ASCII gdal_translate mosaic
TIFF stores floating point values as IEEE754 floats. Talking about significant figures doesn't make much sense. You could test using Float64 with the hope that 50.597 can be exactly represented as a Float64. Otherwise you'll have to do the rounding when reading back from TIFF. Hi Even, Thanks for the information. I don't pretend to understand the fine details of floats, but it does seem counter intuitive that a simple 2 or 3 decimal places can't be accurately represented. I tried the Float64 but the values are identical (16 significant figures) even though the filesize is predictably larger. I guess I'll have to make do, but it does introduce the problem of False Precision. Thanks all, Jonathan -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Extra Decimal places added to ASCII gdal_translate mosaic
Hi Even, Thanks for the information. I tried the Float64 but the values are identical (16 significant figures) even though the filesize is predictably larger. I guess I'll have to make do, but it does introduce the problem of False Precision. Very few formats take into account precision. And the GDAL/OGR model has currently no explicit support for that. My comment about False Precision was referring to the engineering/science problem: https://en.wikipedia.org/wiki/False_precision In our case this is height data, so now folks will see heights down to... femtometres?! You and I know that's obviously wrong, but not all users are quite so technical. Thanks again, Jonathan -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] Gdal GeoTIFFs in MapInfo
Hi List, I've created a GeoTIFF with GDAL. It works fine in QGIS and ArcMap, but MapInfo refuses to load it, giving me a Image file open error message. I've come across several problems and was wondering if anyone else had experience with GDAL created files in MapInfo. Problems: 1) MapInfo can't handle compressed overviews ( *--config COMPRESS_OVERVIEW DEFLATE*) - I can easily get past this one by not compressing them. 2) However, my next problem is that MapInfo can't read the projection information from within the file. Is GDAL capable of creating a TAB file? 3) Finally, while MapInfo will happily open my small subset sample (above issues not withstanding), it won't load my full image (about 2GB). Both are created using the exact same commands: gdalbuildvrt -srcnodata - -vrtnodata - -input_file_list asc_list.txt merged_vrt.vrt gdal_translate -of GTiff -co TILED=YES -co COMPRESS=DEFLATE -co BLOCKXSIZE=512 -co BLOCKYSIZE=512 -a_srs EPSG:27700 merged_vrt.vrt DSM_2012-13.tif gdaladdo -r average DSM_2012-13.tif 2 4 8 16 32 64 It's not the overviews because it still happens when they're not generated. Thoughts welcome, Jonathan -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Extra Decimal places added to ASCII gdal_translate mosaic
Hi List, Bringing back an old thread but with a slight change. I'm now trying to get GDAL to mosaic and export from the ASCII files to a GeoTIFF. But the GeoTIFF writer in GDAL is doing the exact same thing as the ASCII writer was - it's taking *50.597* (source value) and converting it into something like *50.596895342423*. I've tried using both gdal_merge and going via a VRT (gdalbuildvrt - gdal_translate). The result is the same. I also tried the *DECIMAL_PRECISION* trick, but gdal informs me: Warning 6: Driver GTiff does not support DECIMAL_PRECISION creation option How do I get GDAL to merge these into a GeoTIFF without altering the data? Thanks, Jonathan On 24 January 2014 16:12, Jonathan Moules jonathanmou...@warwickshire.gov.uk wrote: Thanks Norman. I'd already seen that page but failed to read it! However, having tested it, it's not functioning as the manual says. As best I can tell, DECIMAL_PRECISION=3 should give me 3 decimal places. But what it instead does it give me 3 significant figures. So 98.354 has become 98.3. This becomes a problem because many of my heights are hovering around 100. So some will either have to be 1 character too long to ensure no data loss for the rest. Is this a bug? The documentation on the page linked to is quite clear it should be decimal places. Less of an issue, but it also didn't make any difference to: xllcorner432000. yllcorner242000. cellsize 2. GDAL 1.10.1, released 2013/08/26 Cheers, Jonathan On 24 January 2014 15:11, Norman Vine n...@cape.com wrote: try gdal_translate -of AAIGrid -co DECIMAL_PRECISION=3 abc.vrt abc.asc see http://www.gdal.org/frmt_various.html HTH Norman On Jan 24, 2014, at 10:01 AM, Jonathan Moules jonathanmou...@warwickshire.gov.uk wrote: Hi List, I'm trying to mosaic some ASCII grid tiles into a single large ASCII file (which can then easily be handled). I'm using this process: REM Create list of files dir /b /s *.asc asc_list.txt REM Turn into VRT gdalbuildvrt -srcnodata - -vrtnodata 0 -a_srs EPSG:27700 -input_file_list asc_list.txt abc.vrt REM actual mosaicing gdal_translate -of AAIGrid abc.vrt abc.asc The input files are like: ncols 500 nrows 500 xllcorner 432000.000 yllcorner 243000.000 cellsize 2 nodata_value -.0 98.354 98.449 98.658 98.874 99.038 99.096 But the outputs are like: fncols500 nrows1000 xllcorner432000. yllcorner242000. cellsize 2. NODATA_value 0 98.353996276855469 98.448997497558594 98.657997131347656 98.874000549316406 Suddenly everything has a dozen extra decimal places. Not only does this induce the False Precision problem (we don't have heights to the femtometer!), it also massively increases the filesize. How can I get this to not happen? Thanks, Jonathan This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Extra Decimal places added to ASCII gdal_translate mosaic
Hi Etienne, I was looking for something like that but the http://www.gdal.org/frmt_gtiff.html page didn't list it explicitly. That said, looking at the gdalinfo now (a little late) I see that it's already defaulting to Float32. I also tested with the -ot Float32 arguement, but there was no change. I'm counting 16 significant figures; I only want 5 or 6. Cheers, Jonathan On 10 March 2014 17:49, Etienne Tourigny etourigny@gmail.com wrote: Have you tried using a 36bit (float) instead of 64bit (double) output? double is usually the default output type and has more significant digits (~16) than your source, so it's representation is ok in double format see [1] for an explanation on floating-point precision the following argument would do this -ot Float32 [1] http://en.wikipedia.org/wiki/Floating_point Etienne On Mon, Mar 10, 2014 at 2:25 PM, Jonathan Moules jonathanmou...@warwickshire.gov.uk wrote: Hi List, Bringing back an old thread but with a slight change. I'm now trying to get GDAL to mosaic and export from the ASCII files to a GeoTIFF. But the GeoTIFF writer in GDAL is doing the exact same thing as the ASCII writer was - it's taking *50.597* (source value) and converting it into something like *50.596895342423*. I've tried using both gdal_merge and going via a VRT (gdalbuildvrt - gdal_translate). The result is the same. I also tried the *DECIMAL_PRECISION* trick, but gdal informs me: Warning 6: Driver GTiff does not support DECIMAL_PRECISION creation option How do I get GDAL to merge these into a GeoTIFF without altering the data? Thanks, Jonathan On 24 January 2014 16:12, Jonathan Moules jonathanmou...@warwickshire.gov.uk wrote: Thanks Norman. I'd already seen that page but failed to read it! However, having tested it, it's not functioning as the manual says. As best I can tell, DECIMAL_PRECISION=3 should give me 3 decimal places. But what it instead does it give me 3 significant figures. So 98.354 has become 98.3. This becomes a problem because many of my heights are hovering around 100. So some will either have to be 1 character too long to ensure no data loss for the rest. Is this a bug? The documentation on the page linked to is quite clear it should be decimal places. Less of an issue, but it also didn't make any difference to: xllcorner432000. yllcorner242000. cellsize 2. GDAL 1.10.1, released 2013/08/26 Cheers, Jonathan On 24 January 2014 15:11, Norman Vine n...@cape.com wrote: try gdal_translate -of AAIGrid -co DECIMAL_PRECISION=3 abc.vrt abc.asc see http://www.gdal.org/frmt_various.html HTH Norman On Jan 24, 2014, at 10:01 AM, Jonathan Moules jonathanmou...@warwickshire.gov.uk wrote: Hi List, I'm trying to mosaic some ASCII grid tiles into a single large ASCII file (which can then easily be handled). I'm using this process: REM Create list of files dir /b /s *.asc asc_list.txt REM Turn into VRT gdalbuildvrt -srcnodata - -vrtnodata 0 -a_srs EPSG:27700 -input_file_list asc_list.txt abc.vrt REM actual mosaicing gdal_translate -of AAIGrid abc.vrt abc.asc The input files are like: ncols 500 nrows 500 xllcorner 432000.000 yllcorner 243000.000 cellsize 2 nodata_value -.0 98.354 98.449 98.658 98.874 99.038 99.096 But the outputs are like: fncols500 nrows1000 xllcorner432000. yllcorner242000. cellsize 2. NODATA_value 0 98.353996276855469 98.448997497558594 98.657997131347656 98.874000549316406 Suddenly everything has a dozen extra decimal places. Not only does this induce the False Precision problem (we don't have heights to the femtometer!), it also massively increases the filesize. How can I get this to not happen? Thanks, Jonathan This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received
[gdal-dev] Extra Decimal places added to ASCII gdal_translate mosaic
Hi List, I'm trying to mosaic some ASCII grid tiles into a single large ASCII file (which can then easily be handled). I'm using this process: REM Create list of files dir /b /s *.asc asc_list.txt REM Turn into VRT gdalbuildvrt -srcnodata - -vrtnodata 0 -a_srs EPSG:27700 -input_file_list asc_list.txt abc.vrt REM actual mosaicing gdal_translate -of AAIGrid abc.vrt abc.asc The input files are like: ncols 500 nrows 500 xllcorner 432000.000 yllcorner 243000.000 cellsize 2 nodata_value -.0 98.354 98.449 98.658 98.874 99.038 99.096 But the outputs are like: fncols500 nrows1000 xllcorner432000. yllcorner242000. cellsize 2. NODATA_value 0 98.353996276855469 98.448997497558594 98.657997131347656 98.874000549316406 Suddenly everything has a dozen extra decimal places. Not only does this induce the False Precision problem (we don't have heights to the femtometer!), it also massively increases the filesize. How can I get this to not happen? Thanks, Jonathan -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Extra Decimal places added to ASCII gdal_translate mosaic
Thanks Norman. I'd already seen that page but failed to read it! However, having tested it, it's not functioning as the manual says. As best I can tell, DECIMAL_PRECISION=3 should give me 3 decimal places. But what it instead does it give me 3 significant figures. So 98.354 has become 98.3. This becomes a problem because many of my heights are hovering around 100. So some will either have to be 1 character too long to ensure no data loss for the rest. Is this a bug? The documentation on the page linked to is quite clear it should be decimal places. Less of an issue, but it also didn't make any difference to: xllcorner432000. yllcorner242000. cellsize 2. GDAL 1.10.1, released 2013/08/26 Cheers, Jonathan On 24 January 2014 15:11, Norman Vine n...@cape.com wrote: try gdal_translate -of AAIGrid -co DECIMAL_PRECISION=3 abc.vrt abc.asc see http://www.gdal.org/frmt_various.html HTH Norman On Jan 24, 2014, at 10:01 AM, Jonathan Moules jonathanmou...@warwickshire.gov.uk wrote: Hi List, I'm trying to mosaic some ASCII grid tiles into a single large ASCII file (which can then easily be handled). I'm using this process: REM Create list of files dir /b /s *.asc asc_list.txt REM Turn into VRT gdalbuildvrt -srcnodata - -vrtnodata 0 -a_srs EPSG:27700 -input_file_list asc_list.txt abc.vrt REM actual mosaicing gdal_translate -of AAIGrid abc.vrt abc.asc The input files are like: ncols 500 nrows 500 xllcorner 432000.000 yllcorner 243000.000 cellsize 2 nodata_value -.0 98.354 98.449 98.658 98.874 99.038 99.096 But the outputs are like: fncols500 nrows1000 xllcorner432000. yllcorner242000. cellsize 2. NODATA_value 0 98.353996276855469 98.448997497558594 98.657997131347656 98.874000549316406 Suddenly everything has a dozen extra decimal places. Not only does this induce the False Precision problem (we don't have heights to the femtometer!), it also massively increases the filesize. How can I get this to not happen? Thanks, Jonathan This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Setting the noData flag in a GeoTIFF
Hi Jukka, Good point. In fact, that's one of the reasons I want to set the noData flag. I figure when compressing hopefully the compressor will see that it is noData and will treat it differently from the way it would an actual value (which results in things like your mentioned 255,255,254) so I won't get the white-border-effect. I guess I'll find out if/when I get to that point. Cheers, Jonathan On 4 December 2013 06:58, Jukka Rahkonen jukka.rahko...@mmmtike.fi wrote: Jonathan Moules jonathanmoules at warwickshire.gov.uk writes: Hi Even, That worked great. It figures that I'd tried about 6 different permutations, including a few without -hidenodata; but none of them was the right one it seems. It's no problem recreating the tifs as they're small test tifs that takes a few seconds to process. Many thanks!Jonathan Hi, Be warned that the result may not be what you believe. You are compressing to jpeg and as a lossy method it does not store constant pixel values for the whole nodata area. If 255-255-255 pixels are nodata values you will get also pixels like 255-255-254 here and there and they will not be transparent in applications. That may not be a problem for you if you combine rectangular images with full of data into a mosaic but if you warp images one by one into another projection so that they will get rotated and then want to use the warped images side by side you would see that the overlapping areas look bad. -Jukka Rahkonen- ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] -wm (memory in mb with gdalwarp)
Thanks for the replies. I think in this instance that as Even suggests it's because both the ECW SDK and GDAL are trying to allocate which will take me over the RAM limit. From a users perspective I just figured that more RAM = better, based on experience that's usually how things work. Useful to know it's not here. Thanks again, Jonathan On 2 December 2013 19:03, Even Rouault even.roua...@mines-paris.org wrote: Le lundi 02 décembre 2013 19:42:59, Jonathan Moules a écrit : On 2 December 2013 18:32, Even Rouault even.roua...@mines-paris.org wrote: Le lundi 02 décembre 2013 19:12:07, Jonathan Moules a écrit : Hi Even, To hit the 10GB limit I'm using (Windows Server 2008): gdal_translate -of GTiff -co TILED=YES -co BIGTIFF=YES -co COMPRESS=JPEG -co JPEG_QUALITY=60 -co BLOCKXSIZE=512 -co BLOCKYSIZE=512 -co PHOTOMETRIC=YCBCR 2000.vrt 2000.tif I tried adding: --config GDAL_CACHEMAX 270 but then it'd just crash with out-of-memory (it wouldn't actually stop at 27GB; so I guess that works differently again). Crash (Windows dialog) or clean exit with a GDAL error message ? What are the files (TIFF, ECW, ...) in the VRT ? I don't see anything wrong, and have not OS and hardware specs to reproduce that. A clean GDAL error: 0...ERROR 2: GDALRasterBlock::Internalize : Out of memory allocating 262144 It is really the operating system refusing to do the allocation. Perhaps there are too much of small allocation, but the virtual memory space on 64 bit is huge so fragmentation should not be an issue. Or perhaps you've just reached the RAM size due to the ECW SDK eating RAM by itself. According to http://www.gdal.org/frmt_ecw.html, it can consume up to 1/4 of the RAM by default. I'm not sure which version of the SDK you are using, but I remember that with the ECW SDK 3.3 it requires a patch for Linux 64 bit, otherwise it could use much more RAM than 1/4. You could try setting ECW_CACHE_MAXMEM explicitely (possibly below 2 GB to avoid any potential 32/64 bit issue) bytes. ERROR 1: GetBlockRef failed at X block offset 334, Y block offset 60 The files are one single ECW file. 5GB compressed, about 200GB uncompressed. -- Geospatial professional services http://even.rouault.free.fr/services.html -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] Setting the noData flag in a GeoTIFF
Hi List, I have a question about setting a value as noData in GeoTIFFs. I'm creating a Geotiff with this two stage process (mosaicing a lot of tiles together): gdalbuildvrt -hidenodata -srcnodata 255 -vrtnodata 255 -input_file_list tiff_list.txt 255.vrt gdal_translate -of GTiff -co TILED=YES -co BIGTIFF=YES -co COMPRESS=JPEG -co JPEG_QUALITY=80 -co BLOCKXSIZE=512 -co BLOCKYSIZE=512 -co PHOTOMETRIC=YCBCR 255.vrt 255mask.tif The resultant image is fine, but when I access it, both via ArcGIS and QGIS it says that there is no noData value set. I've tried reading into it and came across this - http://www.gdal.org/frmt_gtiff.html - but it requires a deeper level of the GeoTIFF library than I have to make much sense of it. I did see GDAL_TIFF_INTERNAL_MASK and tried to set it to true (using -co and --config), but GDAL just complained about it not being a valid configuration option. Can GeoTiffs have a noData value in them that applications will see and honour? If so, can I set it in GDAL? It's not critically important, more nice-to-have and also trying to understand the format better Thanks, Jonathan -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] GDAL python tools and batch files
Hi List, I'm trying to batch (Windows) a fairly simple task - mosaic and then pyramid some rasters. While testing, I've simplified my batch file down to two lines: gdal_merge -o %THIS_DIR%.tif -of GTiff -co TILED=YES -co BIGTIFF=YES -co COMPRESS=JPEG -co JPEG_QUALITY=80 -co BLOCKXSIZE=512 -co BLOCKYSIZE=512 --optfile tiff_list.txt -co PHOTOMETRIC=YCBCR pause The problem is that once gdal_merge has finished processing, rather than continue onward the batch file ends, so in this case the pause is never triggered. Doesn't matter what comes next, the batch ends as soon as the merge is complete. However the same simple batch file but with a gdaladdo command instead works fine. Is this a python thing? Anyone know why that is? Should it be doing this? Using the OSGeo GDAL 1.10.1 build, but it happened with FWTools too (thought it was a FWTools bug at the time). Tested on Windows 7 and Server 2008. Cheers, Jonathan -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] GDAL python tools and batch files
Thanks Paul, that did it. I guess it's a python/batch file interaction thing. Cheers, Jonathan On 2 December 2013 11:51, Paul Meems p.me...@topx-group.nl wrote: I think you should add 'call' at the begining of the command: http://www.microsoft.com/resources/documentation/windows/xp/all/proddocs/en-us/call.mspx?mfr=true Met vriendelijke groet, Paul Meems *Paul Meems*Senior GIS consultant 06-53989481 TopX Geo-ICT http://www.topx-geo-ict.nl http://topx-group.nl/topx-geo-ictWij bieden ondersteuning voor MapWindow GIS http://www.mapwindow.org/ 2013/12/2 Jonathan Moules jonathanmou...@warwickshire.gov.uk Hi List, I'm trying to batch (Windows) a fairly simple task - mosaic and then pyramid some rasters. While testing, I've simplified my batch file down to two lines: gdal_merge -o %THIS_DIR%.tif -of GTiff -co TILED=YES -co BIGTIFF=YES -co COMPRESS=JPEG -co JPEG_QUALITY=80 -co BLOCKXSIZE=512 -co BLOCKYSIZE=512 --optfile tiff_list.txt -co PHOTOMETRIC=YCBCR pause The problem is that once gdal_merge has finished processing, rather than continue onward the batch file ends, so in this case the pause is never triggered. Doesn't matter what comes next, the batch ends as soon as the merge is complete. However the same simple batch file but with a gdaladdo command instead works fine. Is this a python thing? Anyone know why that is? Should it be doing this? Using the OSGeo GDAL 1.10.1 build, but it happened with FWTools too (thought it was a FWTools bug at the time). Tested on Windows 7 and Server 2008. Cheers, Jonathan This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] -wm (memory in mb with gdalwarp)
Hi List, I've just tried the following: call gdalwarp -of GTiff *-wm 13000* -co TILED=YES -co BIGTIFF=YES -co COMPRESS=JPEG -co JPEG_QUALITY=80 -co BLOCKXSIZE=512 -co BLOCKYSIZE=512 -co PHOTOMETRIC=YCBCR SP1937.tif SP1938.tif SP1939.tif 1.tif The documentation says: -wm memory_in_mb: Set the amount of memory (in megabytes) that the warp API is allowed to use for caching. So by my reading I'm allocating a bit less than 13GB (of my 16,730,672kB RAM) to gdalwarp. Imagine my surprise when I get this error: ERROR 5: GDALWarpOptions.Validate() dfWarpMemoryLimit=13000 is unreasonably small. After more testing, anything with 4 digits works fine (including ), but anything with five digits (i.e. 10,000) gives me that error. Is this a bug? (GDAL 1.10.1) Cheers, Jonathan -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] -wm (memory in mb with gdalwarp)
Hi Even, A further thought while we're on the issue of limits, it appears that gdal (both gdal_translate and gdalwarp) have an internal limit of 10GB when they're auto-allocating RAM. Assuming this is the case and it's not just coincidence, it might it be worth updating them to allow more use if appropriate (it seems to be hitting this for a large-ecw to geotiff convesion)? For instance my work desktop has 16GB and the server I'm also doing some processing on has 36GB. Cheers, Jonathan (also I failed to reply-all to my previous reply; so included below) On 2 December 2013 17:25, Jonathan Moules jonathanmou...@warwickshire.gov.uk wrote: Hi Even, Thanks for the detailed reply. It seems in the end that I can use a combination of gdalbuildvrt and gdal_translate to get the same results without having to bother with manual memory management. The output is the same size and the quality is identical. From a user perspective, I would suggest removing this behaviour (though documenting would work too). You could add a gdal warning if they set it above 10GB with a similar caveat to the Note in your email. Of course, then in another 7 years you'll be asked why the warning is triggering at the low 10GB! Thanks again, Jonathan On 2 December 2013 16:39, Even Rouault even.roua...@mines-paris.orgwrote: Le lundi 02 décembre 2013 17:17:11, Jonathan Moules a écrit : Hi List, I've just tried the following: call gdalwarp -of GTiff *-wm 13000* -co TILED=YES -co BIGTIFF=YES -co COMPRESS=JPEG -co JPEG_QUALITY=80 -co BLOCKXSIZE=512 -co BLOCKYSIZE=512 -co PHOTOMETRIC=YCBCR SP1937.tif SP1938.tif SP1939.tif 1.tif The documentation says: -wm memory_in_mb: Set the amount of memory (in megabytes) that the warp API is allowed to use for caching. So by my reading I'm allocating a bit less than 13GB (of my 16,730,672kB RAM) to gdalwarp. Imagine my surprise when I get this error: ERROR 5: GDALWarpOptions.Validate() dfWarpMemoryLimit=13000 is unreasonably small. After more testing, anything with 4 digits works fine (including ), but anything with five digits (i.e. 10,000) gives me that error. Is this a bug? Rather a undocumented attempt to correct wrong user supplied parameter. In gdalwarp utility, if the value of -wm is 1, then it is considered as megabytes ( and multiplied by 1024*1024 for GDAL internals), otherwise it is considered as bytes, and passed unmodified to GDAL internals. (the last change to that logic was in http://trac.osgeo.org/gdal/changeset/10817 , 7 years ago, where the threshold was modified from 4000 to 1, so at that time 10 GB was really enormous !) In gdal warper code, there is a test that checks if the memory limit variable (that must be in bytes now) is at least 10, and bail out if it is not the case. So currently you could specify the size as bytes : -wm 16730672000. I'm wondering if we should not remove this GDAL-isense logic in gdalwarp utility and just implement documented behaviour... Or maybe document the current behaviour. Note: unless your warping process implies considerable image shape distortion (which is unlikely to be the case here since I don't see any -s_srs or -t_srs flag), you don't need such a huge value. It could decrease performance actually. The only case where it might help is if your source images are not tiled. (GDAL 1.10.1) Cheers, Jonathan -- Geospatial professional services http://even.rouault.free.fr/services.html -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] -wm (memory in mb with gdalwarp)
Hi Even, To hit the 10GB limit I'm using (Windows Server 2008): gdal_translate -of GTiff -co TILED=YES -co BIGTIFF=YES -co COMPRESS=JPEG -co JPEG_QUALITY=60 -co BLOCKXSIZE=512 -co BLOCKYSIZE=512 -co PHOTOMETRIC=YCBCR 2000.vrt 2000.tif I tried adding: --config GDAL_CACHEMAX 270 but then it'd just crash with out-of-memory (it wouldn't actually stop at 27GB; so I guess that works differently again). I was also using this (which is a simple tif to tif) and it too used 10GB before stopping itself (this was Windows 7): call gdal_merge -o output.tif -of GTiff -co TILED=YES -co BIGTIFF=YES -co COMPRESS=JPEG -co JPEG_QUALITY=80 -co BLOCKXSIZE=512 -co BLOCKYSIZE=512 --optfile tiff_list.txt -co PHOTOMETRIC=YCBCR It may just have been co-incident. I think I'll just let it do it automatically and hope for the best. Jonathan On 2 December 2013 17:57, Even Rouault even.roua...@mines-paris.org wrote: Le lundi 02 décembre 2013 18:49:42, Jonathan Moules a écrit : Hi Even, A further thought while we're on the issue of limits, it appears that gdal (both gdal_translate and gdalwarp) have an internal limit of 10GB when they're auto-allocating RAM. I'm not aware of such a limitation in the code. It might be that your OS cannot allocate that much memory in one single malloc() if it is the allocation pattern involved here, but I'm just guessing. What commands did you exactly tried to trigger this ? Assuming this is the case and it's not just coincidence, it might it be worth updating them to allow more use if appropriate (it seems to be hitting this for a large-ecw to geotiff convesion)? The ECW SDK itself has its own RAM management. See http://www.gdal.org/frmt_ecw.html For instance my work desktop has 16GB and the server I'm also doing some processing on has 36GB. Cheers, Jonathan (also I failed to reply-all to my previous reply; so included below) On 2 December 2013 17:25, Jonathan Moules jonathanmou...@warwickshire.gov.uk wrote: Hi Even, Thanks for the detailed reply. It seems in the end that I can use a combination of gdalbuildvrt and gdal_translate to get the same results without having to bother with manual memory management. The output is the same size and the quality is identical. From a user perspective, I would suggest removing this behaviour (though documenting would work too). You could add a gdal warning if they set it above 10GB with a similar caveat to the Note in your email. Of course, then in another 7 years you'll be asked why the warning is triggering at the low 10GB! Thanks again, Jonathan On 2 December 2013 16:39, Even Rouault even.rouault@mines- paris.orgwrote: Le lundi 02 décembre 2013 17:17:11, Jonathan Moules a écrit : Hi List, I've just tried the following: call gdalwarp -of GTiff *-wm 13000* -co TILED=YES -co BIGTIFF=YES -co COMPRESS=JPEG -co JPEG_QUALITY=80 -co BLOCKXSIZE=512 -co BLOCKYSIZE=512 -co PHOTOMETRIC=YCBCR SP1937.tif SP1938.tif SP1939.tif 1.tif The documentation says: -wm memory_in_mb: Set the amount of memory (in megabytes) that the warp API is allowed to use for caching. So by my reading I'm allocating a bit less than 13GB (of my 16,730,672kB RAM) to gdalwarp. Imagine my surprise when I get this error: ERROR 5: GDALWarpOptions.Validate() dfWarpMemoryLimit=13000 is unreasonably small. After more testing, anything with 4 digits works fine (including ), but anything with five digits (i.e. 10,000) gives me that error. Is this a bug? Rather a undocumented attempt to correct wrong user supplied parameter. In gdalwarp utility, if the value of -wm is 1, then it is considered as megabytes ( and multiplied by 1024*1024 for GDAL internals), otherwise it is considered as bytes, and passed unmodified to GDAL internals. (the last change to that logic was in http://trac.osgeo.org/gdal/changeset/10817 , 7 years ago, where the threshold was modified from 4000 to 1, so at that time 10 GB was really enormous !) In gdal warper code, there is a test that checks if the memory limit variable (that must be in bytes now) is at least 10, and bail out if it is not the case. So currently you could specify the size as bytes : -wm 16730672000. I'm wondering if we should not remove this GDAL-isense logic in gdalwarp utility and just implement documented behaviour... Or maybe document the current behaviour. Note: unless your warping process implies considerable image shape distortion (which is unlikely to be the case here since I don't see any -s_srs or -t_srs flag), you don't need such a huge value. It could decrease performance actually. The only case where it might help is if your source images are not tiled. (GDAL 1.10.1
Re: [gdal-dev] -wm (memory in mb with gdalwarp)
On 2 December 2013 18:32, Even Rouault even.roua...@mines-paris.org wrote: Le lundi 02 décembre 2013 19:12:07, Jonathan Moules a écrit : Hi Even, To hit the 10GB limit I'm using (Windows Server 2008): gdal_translate -of GTiff -co TILED=YES -co BIGTIFF=YES -co COMPRESS=JPEG -co JPEG_QUALITY=60 -co BLOCKXSIZE=512 -co BLOCKYSIZE=512 -co PHOTOMETRIC=YCBCR 2000.vrt 2000.tif I tried adding: --config GDAL_CACHEMAX 270 but then it'd just crash with out-of-memory (it wouldn't actually stop at 27GB; so I guess that works differently again). Crash (Windows dialog) or clean exit with a GDAL error message ? What are the files (TIFF, ECW, ...) in the VRT ? I don't see anything wrong, and have not OS and hardware specs to reproduce that. A clean GDAL error: 0...ERROR 2: GDALRasterBlock::Internalize : Out of memory allocating 262144 bytes. ERROR 1: GetBlockRef failed at X block offset 334, Y block offset 60 The files are one single ECW file. 5GB compressed, about 200GB uncompressed. -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] How to get GDAL
- From there I've tried both the msi (plus python bindings extra) and the zip Compiled binaries in a single .zip package- both work for the pure-GDAL stuff (gdalinfo) but none of the python scripts works (i.e. gdal_merge.py). I prefer the zip version, expand it to a folder of my choice, and run the SDKshell.bat. In the command window I can do the things I want to do in GDAL, without any interfering of other programms, even not from QGIS. Inside the shell, all varaibles are set correctly, and therefore all python modules run too. And it is no problem installing an experimental new GDAL version and keeping the old one in parallel. That would be my preferred option, but as noted, the python stuff doesn't work. I was getting a python error whenever I tried to run them. Fortunately the OSGeo4W version works fine. Cheers, Jonathan -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] How to get GDAL
Hi List, I'm bringing this one over from the GeoServer-Users list. Basically I can't for the life of me get a good install of GDAL for Windows where everything *just works*. To date I've tried: http://gisinternals.com/SDK/P - But there are far too many options and I can never be sure I'm getting the right thing (the page assumes a significant level of knowledge about what package is what). In the end I just go with - http://gisinternals.com/SDK/PackageList.aspx?file=release-1600-x64-gdal-1-10-mapserver-6-4.zip - From there I've tried both the msi (plus python bindings extra) and the zip Compiled binaries in a single .zip package- both work for the pure-GDAL stuff (gdalinfo) but none of the python scripts works (i.e. gdal_merge.py). I found FWTools - http://fwtools.maptools.org/ - except that it seems to be an outdated version of GDAL so I'm hitting some GDAL bugs there that I think are since fixed. The QGIS version comes with all manner of environmental variable issues that stop it from working stand-alone for me. I think I may have tried the version that comes with FME too, but similar issues. So in short, where can I go to find a pre-compiled (Windows, ideally 64bit) version of GDAL that Just Works? Thanks, Jonathan -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] How to get GDAL
Wow, that actually worked! /me is pleasantly surprised. :-) Thanks! On 27 November 2013 16:15, Mateusz Loskot mate...@loskot.net wrote: On 27 November 2013 16:02, Jonathan Moules jonathanmou...@warwickshire.gov.uk wrote: So in short, where can I go to find a pre-compiled (Windows, ideally 64bit) version of GDAL that Just Works? http://trac.osgeo.org/osgeo4w/ Best regards, -- Mateusz Loskot, http://mateusz.loskot.net -- This transmission is intended for the named addressee(s) only and may contain sensitive or protectively marked material up to RESTRICTED and should be handled accordingly. Unless you are the named addressee (or authorised to receive it for the addressee) you may not copy or use it, or disclose it to anyone else. If you have received this transmission in error please notify the sender immediately. All email traffic sent to or from us, including without limitation all GCSX traffic, may be subject to recording and/or monitoring in accordance with relevant legislation. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
