"and the Anaconda Navigator does not give me higher options" just run the command "conda install -c conda-forge gdal".
It has to come from conda-forge and not conda-main. You can add conda-forge to the navigator - but it is usually just quicker to use the cli. Sorry - but I don't really have the time to debug test files. On Mon, 26 Oct 2020 at 23:00, Nicolas Cadieux <[email protected]> wrote: > Hi, > > I am currently on 3.0.2. and the Anaconda Navigator does not give me > higher options. I will create a new clean environments and try conda > install. I will try changing the code I think my version of the code > probably corresponds to <3 versions of gdal. It's a cut and paste from the > gdal tutorial. Can I send you my test file to test on your code? > > Thanks, I appreciate the help. > > Nicolas > On 2020-10-26 6:31 p.m., Paul Harwood wrote: > > Me too. > > But which conda package for gdal are you using - there are a lot of ones > out there. > > Are you using conda install -c conda-forge gdal ? > > It should be version 3.1.4 at the moment ? > > Also - I don't know if this is part of it but ... > > "srs = osr.SpatialReference() > srs.SetUTM(11, 1) > srs.SetWellKnownGeogCS("NAD27") > dst_ds.SetProjection(srs.ExportToWkt()) > " > > why not just > > "srs = osr.SpatialReference() > srs.SetUTM(11, 1) > srs.SetWellKnownGeogCS("NAD27") > dst_ds.SetSpatialRef(srs) > " > ?? > > On Mon, 26 Oct 2020 at 22:11, Nicolas Cadieux <[email protected]> > wrote: > >> Hi, >> >> Comes back as None for me... :( I'am working with Anaconda. Could it be >> my environment? >> >> Nicolas >> On 2020-10-26 5:04 p.m., Paul Harwood wrote: >> >> GetSpatialRef is OSR based and it is certainly working for me against >> GDAL Raster datasets - it would not be a method against the GDAL Dataset if >> it was only an OGR thing - OGR has its own objects. >> >> This code is working for me at this very minute and giving me the name of >> the CRS for a gtiff >> >> dataset = gdal.Open(str(filepath)) >> if dataset: >> proj = "None" >> crs = dataset.GetSpatialRef() >> if crs: >> proj = crs.GetName() >> return {"type": "gdal", "driver": >> dataset.GetDriver().ShortName, "proj": proj} >> >> On Mon, 26 Oct 2020 at 20:08, Nicolas Cadieux < >> [email protected]> wrote: >> >>> Hi, >>> >>> Thanks for helping. I will put my code at the end. The problem seems >>> to be between gdal python lib 2.2.2 and 3.0.2. Under 2.2.2 and < 3, >>> GetProjection() works well but return nothing with gdal python lib 3.0.2. >>> It's my understanding that .GetSpatialRef is for OGR data (vector). Am I >>> wrong? >>> >>> I was using GetProjection() because that is the method used in the >>> raster_api tutorial. (https://gdal.org/tutorials/raster_api_tut.html) >>> >>> I understand that gdal 3 now relies on proj 6 for the SpatialRef. I'am >>> trying to figure out if is a problem with proj or with gdal... >>> >>> Nicolas >>> >>> >>> code: >>> >>> # -*- coding: utf-8 -*- >>> """ >>> Created on Sun Oct 25 15:09:05 2020 >>> >>> @author: Nicolas >>> """ >>> # >>> https://stackoverflow.com/questions/37648439/simplest-way-to-save-array-into-raster-file-in-python >>> >>> >>> from osgeo import gdal, osr, ogr >>> import numpy >>> >>> print (gdal.__version__) >>> >>> # https://gdal.org/tutorials/raster_api_tut.html >>> fileformat = "GTiff" >>> driver = gdal.GetDriverByName(fileformat) >>> metadata = driver.GetMetadata() >>> if metadata.get(gdal.DCAP_CREATE) == "YES": >>> print("Driver {} supports Create() method.".format(fileformat)) >>> >>> if metadata.get(gdal.DCAP_CREATECOPY) == "YES": >>> print("Driver {} supports CreateCopy() method.".format(fileformat)) >>> dst_ds = driver.Create(r"c:\temp\gdal.tif", xsize=512, ysize=512, >>> bands=1, eType=gdal.GDT_Byte) >>> dst_ds.SetGeoTransform([444720, 30, 0, 3751320, 0, -30]) >>> srs = osr.SpatialReference() >>> srs.SetUTM(11, 1) >>> srs.SetWellKnownGeogCS("NAD27") >>> dst_ds.SetProjection(srs.ExportToWkt()) >>> print ('srs = ',srs)# this is good >>> raster = numpy.zeros((512, 512), dtype=numpy.uint8) >>> dst_ds.GetRasterBand(1).WriteArray(raster) >>> # Once we're done, close properly the dataset >>> srs = None >>> dst_ds = None #srs is file and well georeferenced in Qgis. >>> >>> # >>> --------------------------------------------------------------------------- >>> # open dataset try to read srs >>> # >>> --------------------------------------------------------------------------- >>> raster_ds = gdal.Open(r"C:\temp\gdal.tif") >>> # first try >>> print("Projection is {}".format(raster_ds.GetProjection())) >>> >>> # second try >>> srs = osr.SpatialReference() >>> srs.ImportFromWkt(raster_ds.GetProjectionRef()) >>> print ('srs =', srs) >>> # thrid try >>> srs = osr.SpatialReference() >>> srs = raster_ds.GetProjection() >>> print('srs =', srs) >>> srs = raster_ds.GetProjectionRef() >>> # forth try >>> # from gdal_info >>> https://github.com/OSGeo/gdal/blob/master/gdal/swig/python/samples/gdalinfo.py >>> pszProjection = raster_ds.GetProjectionRef() >>> print(pszProjection) >>> if pszProjection is not None: >>> hSRS = osr.SpatialReference() >>> if hSRS.ImportFromWkt(pszProjection) == gdal.CE_None: >>> pszPrettyWkt = hSRS.ExportToPrettyWkt(False) >>> print("Coordinate System is:\n%s" % pszPrettyWkt) >>> else: >>> print("Coordinate System is `%s'" % pszProjection) >>> >>> >>> On 2020-10-26 3:44 p.m., Paul Harwood wrote: >>> >>> I am not sure why you are using GetProjection and not GetSpatialRef ? >>> Although, the code snippet you included in your email does not seem to use >>> either? >>> >>> I certainly was in the process of working on some code that uses >>> .GeoSptialRef very successfully on GDAL3.1.3 - so I am sure that it is >>> working! >>> >>> On Mon, 26 Oct 2020 at 18:36, Nicolas Cadieux < >>> [email protected]> wrote: >>> >>>> New info! >>>> >>>> GetProjection works with gdal python 2.2.2 and 2.3.3 but not with with >>>> gdal 3.0.2. Is this a change in the gdal library? >>>> >>>> Nicolas >>>> >>>> On 2020-10-26 1:53 p.m., Nicolas Cadieux wrote: >>>> > Hi, >>>> > >>>> > In case mail question was not clear, I have posted the questions >>>> about >>>> > this on stackexchange. I posted a full code that you can run >>>> exposing >>>> > my problem and question. >>>> > >>>> > Thanks >>>> > >>>> > Nicolas >>>> > >>>> > >>>> https://gis.stackexchange.com/questions/377567/cannot-get-projection-from-raster-dataset-using-getprojection >>>> > >>>> > >>>> > On 2020-10-25 9:31 p.m., Nicolas Cadieux wrote: >>>> >> Hi, >>>> >> >>>> >> With the following code, I get an empty string indicating the >>>> >> projection is not valid. >>>> >> >>>> >> from osgeo import gdal, osr >>>> >> raster_ds = gdal.Open(r"C:\temp\180922_WTE3.tif") >>>> >> target_ds = >>>> >> >>>> driver.Create(r"c:\temp\output.tif",xsize=raster_ds_ncol,ysize=raster_ds_nrow,bands >>>> >>>> >> = 1,eType = gdal.GDT_Float32) >>>> >> (array is a numpty array.) >>>> >> raster_srs = osr.SpatialReference() >>>> >> raster_srs.ImportFromWkt(raster_ds.GetProjectionRef()) >>>> >> target_ds.SetProjection(raster_srs.ExportToWkt()) >>>> >> target_ds.GetRasterBand(1).WriteArray(array) >>>> >> raster_ds = None #close dataset >>>> >> target_ds=None >>>> >> >>>> >> Below is the result of gdal info from qgis. File appears to have a >>>> >> valid projection and is properly georeferenced in QGIS using other >>>> >> data sets. Is this projection wrong or am I missing something in >>>> my >>>> >> python code? The goal is to extract the projection from raster_ds >>>> >> set in order to apply to target_ds. I can create the file, apply a >>>> >> geotransform but can't get the projection. Can you help? >>>> >> >>>> >> Thanks/merci! >>>> >> >>>> >> Nicolas >>>> >> >>>> >> QGIS version: 3.14.16-Pi >>>> >> QGIS code revision: df27394552 >>>> >> Qt version: 5.11.2 >>>> >> GDAL version: 3.0.4 >>>> >> GEOS version: 3.8.1-CAPI-1.13.3 >>>> >> PROJ version: Rel. 6.3.2, May 1st, 2020 >>>> >> Processing algorithm… >>>> >> Algorithm 'Raster information' starting… >>>> >> Input parameters: >>>> >> { 'EXTRA' : '', 'INPUT' : 'C:/temp/180922_WTE3.tif', 'MIN_MAX' : >>>> >> False, 'NOGCP' : False, 'NO_METADATA' : False, 'OUTPUT' : >>>> >> 'TEMPORARY_OUTPUT', 'STATS' : False } >>>> >> >>>> >> GDAL command: >>>> >> gdalinfo C:/temp/180922_WTE3.tif >>>> >> GDAL command output: >>>> >> Warning 1: TIFFReadDirectory:Sum of Photometric type-related color >>>> >> channels and ExtraSamples doesn't match SamplesPerPixel. Defining >>>> >> non-color channels as ExtraSamples. >>>> >> >>>> >> Driver: GTiff/GeoTIFF >>>> >> Files: C:/temp/180922_WTE3.tif >>>> >> Size is 1194, 2783 >>>> >> Coordinate System is: >>>> >> PROJCRS["WGS 84 / UTM zone 18N", >>>> >> BASEGEOGCRS["WGS 84", >>>> >> DATUM["World Geodetic System 1984", >>>> >> ELLIPSOID["WGS 84",6378137,298.257223563, >>>> >> LENGTHUNIT["metre",1]]], >>>> >> PRIMEM["Greenwich",0, >>>> >> ANGLEUNIT["degree",0.0174532925199433]], >>>> >> ID["EPSG",4326]], >>>> >> CONVERSION["UTM zone 18N", >>>> >> METHOD["Transverse Mercator", >>>> >> ID["EPSG",9807]], >>>> >> PARAMETER["Latitude of natural origin",0, >>>> >> ANGLEUNIT["degree",0.0174532925199433], >>>> >> ID["EPSG",8801]], >>>> >> PARAMETER["Longitude of natural origin",-75, >>>> >> ANGLEUNIT["degree",0.0174532925199433], >>>> >> ID["EPSG",8802]], >>>> >> PARAMETER["Scale factor at natural origin",0.9996, >>>> >> SCALEUNIT["unity",1], >>>> >> ID["EPSG",8805]], >>>> >> PARAMETER["False easting",500000, >>>> >> LENGTHUNIT["metre",1], >>>> >> ID["EPSG",8806]], >>>> >> PARAMETER["False northing",0, >>>> >> LENGTHUNIT["metre",1], >>>> >> ID["EPSG",8807]]], >>>> >> CS[Cartesian,2], >>>> >> AXIS["(E)",east, >>>> >> ORDER[1], >>>> >> LENGTHUNIT["metre",1]], >>>> >> AXIS["(N)",north, >>>> >> ORDER[2], >>>> >> LENGTHUNIT["metre",1]], >>>> >> USAGE[ >>>> >> SCOPE["unknown"], >>>> >> AREA["World - N hemisphere - 78°W to 72°W - by country"], >>>> >> BBOX[0,-78,84,-72]], >>>> >> ID["EPSG",32618]] >>>> >> Data axis to CRS axis mapping: 1,2 >>>> >> Origin = (459417.200000000011642,5028584.700000000186265) >>>> >> Pixel Size = (0.050000000000000,-0.050000000000000) >>>> >> Metadata: >>>> >> AREA_OR_POINT=Area >>>> >> TIFFTAG_XRESOLUTION=1 >>>> >> TIFFTAG_YRESOLUTION=1 >>>> >> Image Structure Metadata: >>>> >> INTERLEAVE=BAND >>>> >> Corner Coordinates: >>>> >> Upper Left ( 459417.200, 5028584.700) ( 75d31' 7.03"W, 45d24'34.58"N) >>>> >> Lower Left ( 459417.200, 5028445.550) ( 75d31' 6.99"W, 45d24'30.07"N) >>>> >> Upper Right ( 459476.900, 5028584.700) ( 75d31' 4.28"W, >>>> 45d24'34.59"N) >>>> >> Lower Right ( 459476.900, 5028445.550) ( 75d31' 4.24"W, >>>> 45d24'30.08"N) >>>> >> Center ( 459447.050, 5028515.125) ( 75d31' 5.63"W, 45d24'32.33"N) >>>> >> Band 1 Block=1194x1 Type=UInt16, ColorInterp=Red >>>> >> ... >>>> >> Band 288 Block=1194x1 Type=UInt16, ColorInterp=Undefined >>>> >> >>>> >> >>>> >> >>>> _______________________________________________ >>>> gdal-dev mailing list >>>> [email protected] >>>> https://lists.osgeo.org/mailman/listinfo/gdal-dev >>> >>>
_______________________________________________ gdal-dev mailing list [email protected] https://lists.osgeo.org/mailman/listinfo/gdal-dev
