1

I have been working on this iterator that would count instances where data conformed to a set of conditions over a number of rasters. I had people help me with other aspects of this script (Counting events from a GRIB raster), but now I have hit an error where the return is "ValueError: array larger than output file, or offset off edge", Which I can't seem to be able to work arround. The code is provided bellow. It requires the user to input the path to folder and provide a name for the output GeoTIFF raster. Files read are GRIB.

import gdal
import numpy
import osr
import os

del_mapa = str(input('Vpiši pot do datotek v obliki "C:\\Mapa1\\Mapa2\\...\\Zadnja mapa": '))
os.chdir(del_mapa)

pix_vel = 1000
#Defines pixel size
ime_rez = str(input("Poimenuj datoteko z rezultati: ") + ".tif")
#User can provide a name for the output picture
format = "GTiff"
driver = gdal.GetDriverByName(format)

vir_mere = gdal.Open('inca_20140101-0100.grb')
trans = vir_mere.GetGeoTransform()
geotrans0 = trans[0]
geotrans1 = trans[1]
geotrans3 = trans[3]
geotrans5 = trans[5]
#Extracts transformation data from one of the GRIBs
dst_ds = driver.Create(ime_rez, 401, 301, 1, gdal.GDT_CInt32)
dst_ds.SetGeoTransform([geotrans0, geotrans1, 0, geotrans3, 0, geotrans5])
#and pastes it directly into the newly created GeoTIFF
projInfo = vir_mere.GetProjection()
srs = osr.SpatialReference()
srs.ImportFromWkt(vir_mere.GetProjectionRef())

dst_ds.SetProjection(srs.ExportToWkt())

raster = numpy.zeros((401, 301), dtype=numpy.uint8)
dst_ds.GetRasterBand(1).WriteArray(raster)

def krog(baza,novi):    
    x = 0
    while x < 25:
        baza = gdal.Open(ime_rez)
        if x < 10:
            novi = gdal.Open('inca_20140101-%s00.grb') % ('0' + str(x))
            #This should open the next file in line
            if novi > 15 and novi < 25:
                baza = baza + 1
            else:
                baza = baza + 0

        elif x > 10:
            # novi = gdal.Open('inca_20140101-%s00.grb') % (str(x))
            if novi > 15 and novi < 25:
                baza = baza + 1
            else:
                baza = baza + 0
        x += 1 
dst_ds.GetRasterBand(1).WriteArray(baza)        
dst_ds = None
  • 1
    Try: raster = numpy.zeros((301, 401), dtype=numpy.uint8) – xunilk Jul 26 '17 at 14:23
  • Thank you! That fixed the error, but now I can only get a blank raster out of the program and I know that it should produce some hits (there should be more than 0 everywhere) and I tried different conditions... – peroman200 Jul 26 '17 at 17:36
  • What is your dataset raster: raster or baza? I cannot see how baza dataset values were changed for obtaining a new raster. – xunilk Jul 26 '17 at 17:58
  • raster is the one that is filled with numpy zeroes, baza is defined right after the while x < 25: line. – peroman200 Jul 26 '17 at 18:01
  • I'll take a look later. – xunilk Jul 26 '17 at 18:04
0

As @xunilk indicated in a comment. This ValueError is due to a mix up of columns and rows. You were trying to write a 401x301 array into a 301x401 raster.

Your 2nd issue is that you are trying to write an array of numpy.uint8 into a raster of gdal.GDT_CInt32. The resulting image will be no good; often resulting in an all-black raster.

You should use numpy.int32 to match gdal.GDT_CInt32 or else use one of the other gdal data types:

 'GDT_Byte', 'GDT_CFloat32', 'GDT_CFloat64', 'GDT_CInt16', 'GDT_CInt32', 'GDT_Float32', 'GDT_Float64', 'GDT_Int16', 'GDT_Int32', 'GDT_TypeCount', 'GDT_UInt16', 'GDT_UInt32'

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.