CODEWITHSUNDEEP
arraysdimensionalnetcdf4
alice
alice

Reputation: 325

how to count the netCDF for every daily file in python

i want to achieve calculate the number of event over the min of daily netCDF4 file temperature. I have a code like below, but it keep telling me the index out of bounds. The netCDf4 file a three dimensional array with 349x277x2920. the third dimension is time, the tempature was taken every three hours, so total 365 days of temptures. i want to calculate the daily min where the pixel below 266, the count and then count the daily count of all 365 days together.the code is like this:

from osgeo import gdal
from osgeo import osr
import numpy as np
import netCDF4 as net



finalgrid=np.zeros(shape=[277,349],dtype=int)
infile=net.Dataset("W:/air.2m.2014.nc")
value=infile.variables['air'][:,:,]
infile.close()
dailymean=np.full([277,349],300,dtype=float)


for k in range(0,2921,8):
    emptygrid=np.zeros(shape=[277,349],dtype=int)
    for i in range(0,278):
        for j in range(0,350):
            try:
                if value[k][i][j]>1:
                    dailymean[i][j]=min(value[k][i][j],value[k+1][i][j],value[k+2][i][j],value[k+3][i][j],value[k+4][i][j],value[k+5][i][j],value[k+6][i][j],value[k+7][i][j])

            except:
                continue
    emptygrid[dailymean<=266.0]=emptygrid[dailymean<=266.0]+1
    finalgrid=finalgrid+emptygrid
    print k


driver=gdal.GetDriverByName('gtiff')
print finalgrid.max()
outDs=driver.Create('W:/2014aircount.tif',349,277,1,gdal.GDT_Float64)
outBand=outDs.GetRasterBand(1)
outBand.WriteArray(finalgrid,0,0)
ds=gdal.Open("W:/air.2m.2014.nc")
gt = ds.GetGeoTransform()
outDs.SetGeoTransform(gt)
outDs.SetProjection(ds.GetProjection())
#outDs.SetProjection(ds.GetProjection())
outBand.FlushCache()
ds=None

Upvotes: 0

Views: 484

Answers (1)

N1B4
N1B4

Reputation: 3453

There are a few issues with your code, including the read-in of air and slices of air to compute daily min temps in the looping section. Below is an alternative that will also only involve 1 for-loop over the time dimension. 

import netCDF4
import numpy as np

ncfile = netCDF4.Dataset('/nas/home/air.2m.2014.nc', 'r')

# Read in all three dimensions (lat x lon x time) below 
temp = ncfile.variables['air'][:,:,:] 

# Compute the dimension sizes
nlat, nlon, ntim = np.shape(temp)

# Compute the number of days, knowing that the time step is 3-hourly (or 8 snapshots per day)
ndays = int(ntim/8) 

# Create an empty 3D array (time x lat x lon) that will store daily min temps 
daily_min_temp = np.zeros([ndays, nlat, nlon], dtype='f4')

# Compute daily min temps across the grid
for day in range(ndays):
    daily_temp = temp[:,:,day*8:day*8+8]
    daily_min_temp[day,:,:] = np.min(daily_temp, axis=2)

Upvotes: 0

Related Questions