Re: [gdal-dev] ogr2ogr failing on postgis?
Le dimanche 10 avril 2016 20:39:25, Paolo Cavallini a écrit : > Il 10/04/2016 18:01, Even Rouault ha scritto: > > Yes, this is expected and a rather common problem with shapefiles. > > Thanks a lot for the explanation. QGIS therefore creates an incorrect > prj, but I did not noticed it because it loaded it with the correct > EPSG. Ticket upstream opened: > http://hub.qgis.org/issues/14655 No, QGIS did it fine (actually I suspect it generates .prj file through OGR itself). ESRI .prj file are supposed *not* to contain EPSG codes unfortunately. This is the ESRI WKT variant. Which makes it hard to deal with them on the reading side. The issue is more on the OGR side. -- Spatialys - Geospatial professional services http://www.spatialys.com ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] ogr2ogr failing on postgis?
Il 10/04/2016 18:01, Even Rouault ha scritto: > Yes, this is expected and a rather common problem with shapefiles. Thanks a lot for the explanation. QGIS therefore creates an incorrect prj, but I did not noticed it because it loaded it with the correct EPSG. Ticket upstream opened: http://hub.qgis.org/issues/14655 All the best, and thanks. -- Paolo Cavallini - www.faunalia.eu QGIS & PostGIS courses: http://www.faunalia.eu/training.html ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] ogr2ogr failing on postgis?
Hi Paolo, > Hi all, > I have an ogr2ogr command, from shp to pg, that fails because it > requests write access to spatial_ref_sys. An apparently identical > command on another file runs smoothly. Any hint? > Details and data documented here: > http://hub.qgis.org/issues/14650 Yes, this is expected and a rather common problem with shapefiles. As the shapefile .prj doesn't contain any explicit EPSG code, the PostGIS driver tries to find a match in spatial_ref_sys based on an exact match on the WKT, but as the OGC WKT built from the .prj is different from the one in spatial_ref_sys (one of the difference is the absence of the AUTHORITY["EPSG","3003"] node of course !), it must fallback to creating an ad-hoc entry, hence the need for write access. Workaround: provide explicit SRS, in that case -a_srs EPSG:3003 (as I can see in the .qpj) A more "fuzzy" indentification of SRS could help for that. Even -- Spatialys - Geospatial professional services http://www.spatialys.com ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] ogr2ogr failing on postgis?
Hi Paolo, It seems to me that your shape spatial reference not present in spatial_ref_sys and GDAL try to import it into spatial_ref_sys. Your shape file have such SRS: PROJCS["Monte_Mario_Italy_zone_1",GEOGCS["GCS_Monte Mario",DATUM["D_Monte_Mario",SPHEROID["International_1924",6378388,297]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",9],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",150],PARAMETER["false_northing",0],UNIT["Meter",1]] If I execute such SQL: select * from spatial_ref_sys where srtext like '%Monte_Mario%' I received: 4265;"EPSG";4265;"GEOGCS["Monte Mario",DATUM["Monte_Mario",SPHEROID["International 1924",6378388,297,AUTHORITY["EPSG","7022"]],TOWGS84[-104.1,-49.1,-9.9,0.971,-2.917,0.714,-11.68],AUTHORITY["EPSG","6265"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.01745 (...)" 4806;"EPSG";4806;"GEOGCS["Monte Mario (Rome)",DATUM["Monte_Mario_Rome",SPHEROID["International 1924",6378388,297,AUTHORITY["EPSG","7022"]],TOWGS84[-104.1,-49.1,-9.9,0.971,-2.917,0.714,-11.68],AUTHORITY["EPSG","6806"]],PRIMEM["Rome",12.452333,AUTHORITY["EPSG","8906"] (...)" 3003;"EPSG";3003;"PROJCS["Monte Mario / Italy zone 1",GEOGCS["Monte Mario",DATUM["Monte_Mario",SPHEROID["International 1924",6378388,297,AUTHORITY["EPSG","7022"]],TOWGS84[-104.1,-49.1,-9.9,0.971,-2.917,0.714,-11.68],AUTHORITY["EPSG","6265"]],PRIMEM["Greenwich",0,AUTHORITY[" (...)" 3004;"EPSG";3004;"PROJCS["Monte Mario / Italy zone 2",GEOGCS["Monte Mario",DATUM["Monte_Mario",SPHEROID["International 1924",6378388,297,AUTHORITY["EPSG","7022"]],TOWGS84[-104.1,-49.1,-9.9,0.971,-2.917,0.714,-11.68],AUTHORITY["EPSG","6265"]],PRIMEM["Greenwich",0,AUTHORITY[" (...)" 26591;"EPSG";26591;"PROJCS["Monte Mario (Rome) / Italy zone 1 (deprecated)",GEOGCS["Monte Mario (Rome)",DATUM["Monte_Mario_Rome",SPHEROID["International 1924",6378388,297,AUTHORITY["EPSG","7022"]],TOWGS84[-104.1,-49.1,-9.9,0.971,-2.917,0.714,-11.68],AUTHORITY["EPSG","6806"]], (...)" 26592;"EPSG";26592;"PROJCS["Monte Mario (Rome) / Italy zone 2 (deprecated)",GEOGCS["Monte Mario (Rome)",DATUM["Monte_Mario_Rome",SPHEROID["International 1924",6378388,297,AUTHORITY["EPSG","7022"]],TOWGS84[-104.1,-49.1,-9.9,0.971,-2.917,0.714,-11.68],AUTHORITY["EPSG","6806"]], (...)" I think you need to reproject shape file to EPSG: 3003 or defined it for shape file. Best regards, Dmitry 10.04.2016 18:23, Paolo Cavallini пишет: Hi all, I have an ogr2ogr command, from shp to pg, that fails because it requests write access to spatial_ref_sys. An apparently identical command on another file runs smoothly. Any hint? Details and data documented here: http://hub.qgis.org/issues/14650 Thanks. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Map algebra
08.04.2016, 19:28, Gregory, Matthew kirjoitti: I would want any map algebra library to handle three types of issues when dealing with raster data: 1) Handling mask or nodata values (e.g. an analysis mask should be able to be set on any algebra operation that retains or combines nodata values from its component operands) GDAL supports nodata values and mask bands to identify cells that do not contain valid data. In arithmetics I think it is straightforward that x + nodata = nodata. In logical operations it depends on the interpretation of nodata. If it means unknown, then we are dealing with three-valued logics where for example unknown OR true = true. Think about for example a raster which shows polluted areas as a result of two studies (true, false, unknown). If one study shows a site polluted and the other has no data, then the result should be that the site is polluted. I haven't yet got to the point in my code where I look at nodata and mask bands. 2) Handling different spatial extents (e.g. the ability to specify an output extent -- coming from a specific raster or from the union/intersection of rasters) This may be desirable in some cases. For example if we have a huge raster and want to study some parts of it. In general I would, however, leave this for other already existing GDAL tools, especially gdal_translate, whose functionality is available in APIs as well. 3) Handling different cell sizes (e.g. the ability to specify an analysis cell size and resample/reproject(??) all operands based on this parameter) I think using gdal_translate in this case as a first step makes even more sense. Alex asked about how I iterate over blocks when there are two or more bands involved. The basic idea is that all methods operate on one band, and other bands (in my case if there is a second band) are adjusted to that. The algorithm is the following. "band" is a wrapper (struct) for a GDALRasterBand object and block cache. A block cache is simply an array of pointers to data, which are obtained from the GDALRasterBand object with ReadBlock() and which may be written back with WriteBlock(). band1 = gma_band_initialize(b1); band2 = gma_band_initialize(b2); while (iterate) { block_index i; for (i = blocks in band1) { add block i to cache in band1; block1 = get block i from band1; update cache in band1 to allow desired focal distance; update cache in band2 to allow desired focal distance; call method specific callback with (band1, block1, band2, retval, arg); make a note if the callback indicates a need for iteration; based on the return value from the callback return with error, continue with next block from band1, or write block1 to band1 and continue with next block from band1; } } if iteration, then in some cases some band level things need to be done; } empty cache in band1; empty cache in band2; } The cache update function reads all needed blocks into the cache and discards those that are not needed. If there were no discards, the cache would eventually contain the whole band. In the method callback values from the band2 are obtained by first computing the global cell index and then the respective block and cell index within the block. This may not be the best possible algorithm but to me it is pretty understandable, which is a goal in itself. The code is at https://github.com/ajolma/gdal/tree/trunk/gdal/map_algebra The curse of dimensionality comes in at the point, where the above function is called. The method callback is a template function and it needs the datatype of each band that is given to it as an argument. (I see now that I have made also the above function a template function, which may not be needed). Best, Ari ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev