New answers tagged

0

I have solved this issue by changing the driver string which is : RasterProcessing.Convert(@"E:\@@MapInfo_1\Output\2-0.tif", "E:\\@@MapInfo_1\\Jonathan_Test.tif", "MI_GeoTiff_IMG");


0

apologies for all these libraries, i'm not sure exactly which is needed... Good luck! this is kind of complicated. and not exactly sure how it works. library(igraph) library(RColorBrewer) library(SpatialTools) library(sp) library(ANN) library(rgeos) library(maptools) library(sp) library(rgeos) library(plyr) library(FNN) library(tidyr) library(foreign) ...


1

Note the difference between band rendering done by QGIS and the actual raster extent and values. The values of your raster remain the same, regardless of the band rendering (Multicolor band, Singleband gray, Paletted unique/values) you select in Layer properties. If you take your original layer (i.e. HDF4_EOS:EOS_GRID:"MCD19A2.A2018312.h24v06.006....


0

MapInfo .tab files are just text files; it is quite easy to create them yourself. A .tab file for a registered image looks like this: !table !version 300 !charset WindowsLatin1 Definition Table File "test.png" Type "RASTER" (10000,7500) (0,0) Label "top left", (20000,7500) (256,0) Label "top right", (10000,0) (0,192) Label "bottom left" ...


1

You can either use an image layer type and specify a single image with a bounds or you can use a raster layer type where your raster is using the tiling scheme per https://en.wikipedia.org/wiki/Tiled_web_map or Raster images uploaded to mapbox.com as Tilesets. Specifying a bounds on a raster source will simply limit the Mapbox GL JS client to only request ...


3

The first warning means that the value of GeoASCIIParams tag is not read as it was written because the original image is having NULL character in the value of the tag. NULL can be used as a delimiter between strings http://freeimage.sourceforge.net/fnet/html/A633E9A9.htm but obviously GDAL takes just the first string. The second error means that the writer ...


2

You can suppress the warnings (as long as you're sure there's no real issue with your data) with gdal.PushErrorHandler('CPLQuietErrorHandler'). If you do this any errors will also not get printed, so make sure you tell GDAL to raise a Python exception when an error occurs with gdal.UseExceptions(). E.g. # Stop GDAL printing both warnings and errors to ...


0

You can do the conversion from radians to degrees in python. If you combine that with the keyword argument abilities in Proj, then you can do it like so: import math import pyproj myprojection = pyproj.Proj( proj="gnom", lat_0=math.degrees(0.017453292519943295), lon_0=math.degrees(0.017453292519943295), R=1., )


0

Rasters are always rectangular, with their edges oriented with the axes of their coordinate reference system. Georeferencing an image usually rotates that image. If you rotate a rectangular image, you then have to add a nodata border around it to make the output into a rectangle with its edges oriented with its coordinate reference system. The same is true ...


2

One approach is to use the SQLite SQL dialect. It seems to work at least with a point shapefile where I digitized a few points, some of them in the same location. ogrinfo -dialect sqlite -sql "select geometry, count(*) from duplicate_points group by geometry" duplicate_points.shp Layer name: SELECT Geometry: Unknown (any) Feature Count: 4 Extent: (271....


0

My idea is to order point coordinates by its longitude and latitude in this way: Sort by longitude, in for loop iterate over the coordinates while longitude ascends; Then iterate while latitude descends; etc... Then merge lists with coordinates and create polygon. Not pretend this is the most elegant solution, maybe you send your data so I will test my ...


0

I ran into this problem and eventually realized that the gdal version I had in my build.gradle file (v3.0.0) was newer than and apparently incompatible with the gdal version installed on my computer (v2.4.2). I changed the version number in build.gradle to 2.4.2 and it worked.


1

Global Mapper natively reads PDF and geopdf and writes GeoPDF. You can georeference and save in whatever format you want


0

Thanks to user30184 I installed GDAL 2.4.2 instead of 3.x and now the command gdal2tiles run without error.


1

You can convert x,y to col, row coordinates with gdal.InvGeoTransform, and gdal.ApplyGeoTransform. You can reproject coordinates with osr.CoordinateTransformation.TransformPoint. E.g. from osgeo import gdal, osr ds = gdal.Open(someraster) gt = ds.GetGeoTransform() # Geotransforms allow conversion of pixel to map coordinates crs = ds.GetProjection() ...


3

No, it's not a bug. 5487115.2521567 is correct. The GeoTransform doesn't give you the lower left, it gives you the origin. Which is the upper left. You can calculate the lower left from the origin + y pixel size * no. rows: yllcorner = dem_transform [3] + dem_transform[5] * dem_file.RasterYSize print(yllcorner) 5480455.0079735


0

I want to post some additional info for people who might run into problems installing the GDAL extensions on Ubuntu 16.04, following the instructions on https://docs.geoserver.org/stable/en/user/data/raster/gdal.html These are my observations: In my case, setting the GDAL_DATA environment variable was not needed. The gdal data (http://sourceforge.net/...


0

I would expect that in the line: poVRTDS = poDriver->CreateCopy("media/B02_stack.vrt", poSrcDS,FALSE,APPEND_SUBDATASET=YES,NULL,NULL); you need to place quotes around (at least) APPEND_SUBDATASET=YES: poVRTDS = poDriver->CreateCopy("media/B02_stack.vrt", poSrcDS,FALSE,"APPEND_SUBDATASET=YES",NULL,NULL);


3

For the stable series of GeoServer (currently 2.15.x) you can use these instructions: https://docs.geoserver.org/stable/en/user/data/raster/gdal.html The properly tested packages are the ones of the custom GDAL 1.9.2 build, there are hints on how to use a different/newer version too, but it's not recommended.


0

In R you can use rotate function library(raster) library(gdalUtils) workdir <- "Your workind dir" setwd(workdir) ncfname <- "adaptor.mars.internal-1563580591.3629887-31353-13-1b665d79-17ad-44a4-90ec-12c7e371994d.nc" # get the variables you want dname <- c("v10","u10") # open using raster datasetName <-dname[1] r <- raster(ncfname, ...


0

Aren't your path incorrect? I believe you would want something like "gdal_translate -b 1 " raster_input_dir + "/" + raster_filename + " " + outputDir + "/" + lyr.name() + "test.ecw"


0

I tried this string in pyproj: Proj(proj='lcc', R=6371200, lat_1=25, lat_2=25,lon_0=-95, ellps='clrk66') and my results were strange. when I changed the lat field names: Proj(proj='lcc', R=6371200, lat_0=25, lat_1=25,lon_0=-95, ellps='clrk66') I got the results I was expecting.


0

I have encountered the same shift issue. This shift is not limited to GDAL/QGIS. It is also present in ArcMap, MapInfo and Global Mapper (maybe gdal is used in the backend...). A work around is to convert the GeoTiff to an ERS, fix the easting and northing in the .ers file (open with a text editor) and convert the ERS back to a GeoTiff. The conversion ...


0

I also want to comment on how the variables are being read from the netcdf file. read varaiblaes as follows: import gdal import netCDF4 as nc ncfile = nc.Dataset("C:/path_to_file/file.nc", "r") snow = ncfile.variables["snc"][:,:] lat = ncfile.variables["lat"][:] lon = ncfile.variables["lon"][:]


0

This code requires minor changes and it will work. Instead of: dst_ds = gdal.GetDriverByName('GTiff').Create('output.tif', ny, nx, 1, gdal.GDT_Float32) Use: dst_ds = gdal.GetDriverByName('GTiff').Create('output.tif', nx, ny, 1, gdal.GDT_Float32)


2

A couple of things here: Firstly don't set gdata to None until you are COMPLETELY done with that raster, data will return a broken dataset if you do. Please see for a more detailed explanation: https://trac.osgeo.org/gdal/wiki/PythonGotchas Next thing uninstall gdal using pip then reinstall it using the .whl file from here: https://www.lfd.uci.edu/~...


1

With some help from user2856s answer I figured it out. Like he said first I needed to scale down to 0-255 before converting to Byte. Then the issue was I was not exporting the KML along with the tiles... silly mistake. I just put the kml in a folder with all the tiles. Hopefully this answer helps someone else.


0

I managed to figure out what was going on here after several emails with their support team. The projection is indeed a simple geographic lat/long projection with units in decimal degrees. WGS84 is the datum. The extent labels are actually wrong. The units are not in meters they are thankfully in degrees, minutes, seconds. 'UpperLeftPointMtrs=(20000000....


1

Here is a function to convert a QgsRasterLayer to a numpy array without GDAL through using the block method of the QgsRasterDataProvider(link): from numpy import array def convertRasterToNumpyArray(lyr): #Input: QgsRasterLayer values=[] provider= lyr.dataProvider() block = provider.block(1,lyr.extent(),lyr.width(),lyr.height()) for i in ...


0

to output to the output folder, add folder name like this for %f in (*.kml) do ogr2ogr -f "ESRI shapefile" "output\%f.shp" %f Shapefile is a multi-file format, and OGR updates the files in a non-sequential way. You have to zip it afterwards. Merging can also be done using -append,-update. Refer to this


0

Here's another way to do it: import geopandas as gpd import numpy as np # load an example polygons geodataframe gdf_polys = gpd.read_file(gpd.datasets.get_path('nybb')) It looks like the following: # find the bounds of your geodataframe x_min, y_min, x_max, y_max = gdf_polys.total_bounds # set sample size n = 100 # generate random data within the bounds ...


0

In addition to installing GDAL, you also need to install the GDAL Python bindings. You do this with the command pip install gdal I don't know much about your MacOS environment and how you installed Python, so if this does not work, you would need to adjust your environment settings first to make pip work properly. I am also not sure why you are trying to ...


1

Complete workaround for converting each feature to an individual layer and then rasterize each one it can be observed in following code (where my own paths to layers were used): from osgeo import gdal from osgeo import ogr from osgeo import gdalconst import os raster_file = "/home/zeito/pyqgis_data/mask.tif" shape_file = "/home/zeito/pyqgis_data/test2.shp" ...


1

Your code is perfect. I think you are not working in a projection in meters e.g mollweide. All you need to do is reproject the layer to -- +proj=moll +lon_0=0 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m no_defs Then run the code again and it would be 100%


0

It's probably worth adding that if you are applying that SQL to one table in a large database, you can likely speed up the ogr2ogr operation by feeding fake info to the OCI. Let me explain... Here, I direct the OCI to search for a non-existent table called "junk" in DATABASE ogr2ogr -a_srs epsg:3005 -f "KML" X:\test.kml OCI:user/password@DATABASE:junk -...


0

I was also unable to create the polygon fill using the "KML" driver in ogr2ogr, but was able to use a workaround (text-editing the KML afterwards). Process: 1. Create KML via ogr2ogr in the OSGeo4W Shell: ogr2ogr -a_srs epsg:3005 -f "KML" X:\test.kml OCI:user/password@DATABASE:junk -sql @X:\style.sql -dsco NameField="FIRE_NUMBER" -nln "wildfires" -...


2

It might be a problem with your usage of quotes. Try this (single quotes removed). gdal_translate -of GTiff HDF5:"OMPS-NPP_NMTO3-L2_v2.1_2019m0724t000444_o40089_2019m0724t020831.h5"://ScienceData/UVAerosolIndex omps.tif Also make sure your gdal supports the HDF5 format with gdalinfo --formats


0

I resolved this by reinstalling a previous gdal version supported by django (gdal204 instead of gdal300), removing the system enviroment variables refered to gdal and including these in the django project (settings.py file): os.environ['GDAL_DATA'] = r"C:\Python37\Lib\site-packages\osgeo\data\gdal" os.environ['PROJ_LIB'] = r"C:\Python37\Lib\site-packages\...


2

https://trac.osgeo.org/gdal/ticket/6986 provides further details around this issue. It suggests the following as a potential resolution... "This seems specific to jpeg-in-tiff compression, as using COMPRESS=DEFLATE doesn't show the same issue".


1

What I did when I needed to sample the spectral profile of a web map layer was to right click the layer, export the Raw Data to GeoTiff, choose an appropriate extent, and set an appropriate resolution for the layer (it likely will not set a value other than 0, which will yield a failed export). Then, you do you work off this rendered image. Not ideal, but ...


2

I agree with @Devin Cairns. Most likely the situation is a difference between both visualizations rather than their underlying values. If you want to be sure that the values are different inside your area of interest, try the following code: data = gdal.Open(PR_files[10]) img = data.GetRasterBand(1) raster = img.ReadAsArray() nd_value = img.GetNoDataValue() ...


1

I suspect the array you read from the GeoTiff is correct. The reason it appears that you lost all valid data is because of the no data value washing out the valid data in the matplotlib plot colormap. Raster files have an inherent no data or "nulling" value that GIS software implicitly ignores for display and calculations. In your case, -9999 is the no data ...


1

You can use Dataset.RasterCount to test if the dataset is raster and Dataset.GetLayerCount() for vector. Source: http://osgeo-org.1560.x6.nabble.com/Discover-whether-a-GDALDataset-is-raster-or-vector-td5270223.html However, some datasets can contain both raster and vector data... So depending on your usecase you may need to account for this: if ds....


0

I would just match the name of the driver to a simple look up table. Most geospatial applications will only be working with a couple expected data types, so this is pretty easy to implement. from osgeo import gdal dtypes = { 'raster': ['GTiff', 'JPEG'], 'vector': ['ESRI Shapefile', 'GeoJSON'] } def check_type(infile): ds = gdal.OpenEx(infile) ...


1

This might not be the most desired way, but I have found it to be reliable in answering this question. There may be edge cases or drivers that can handle both vector and raster formats. The short answer is to investigate the driver of the opened GDALDataset, and particularly the driver metadata which may include the fields DCAP_RASTER and DCAP_VECTOR. I ...


2

It seems possible that you have an old variable for the array. I cleaned up your code a little by referring directly to arr. ds = gdal.Open(template_raster) band = ds.GetRasterBand(1) arr = band.ReadAsArray() [rows,cols] = arr.shape driver = gdal.GetDriverByName('GTiff') out_ds=driver.Create(output_raster,cols,rows,1,gdal.GDT_Float32) out_band = out_ds....


1

There is a DISCARD_LSB option for GeoTIFFs that could be useful in conjunction with compression. For example the following call will compress input.tif to output.tif using the DEFLATE algorithm with a floating point predictor (2) and discard the 2 least significant bits: gdal_translate input.tif output.tif -co COMPRESS=DEFLATE -co PREDICTOR=2 DISCARD_LSB=2


0

It looks like gdal_rasterize is not in your path, or perhaps GDAL is not installed? I don't have a Mac myself, but it looks like brew install gdal --HEAD (from https://gis.stackexchange.com/a/165286/2397) should install GDAL command line utilities for you if they aren't installed already.


1

out_shape can not have floating point values as numpy arrays can't have floating point shapes. Either explicitly truncate or round the out_shape values or use integer division. Also, note that rasterio reads data in (bands, rows, cols) order and you are requesting a (rows, cols, bands) order array. with rasterio.open('image.tif') as dataset: data = ...


2

dst.write(data) is the way to write a 3D array to raster. However, rasterio is expecting data in (bands, rows, cols) order and you are passing a (rows, cols, bands) order array. Try data = dataset.read( out_shape=(dataset.count, int(dataset.height / 3.75), int(dataset.width / 3.75)), resampling=Resampling.cubic ) Also note that your transform is ...


Top 50 recent answers are included