Latest Contributions by DuminduJayasekera

Peter, I have two questions:  1) I am trying to use the above code with LANDSAT data (LANDSAT 8 OLI_TIRS C1 Level-1) to calculate NDVI for a given state in US. I am getting the following error:  Runtime error Traceback (most recent call last): File "<string>", line 4, in <module> ImportError: No module named pathlib. I tried installing pathlib but i cannot install due to security reasons in my official laptop.i tried os.path but did not work too. Do not know the reason. 2). Is there a work around for the above code for using LANDSAT ARD data to calculate NDVI values.?  Any help is appreciated.  ... View more

Dan and Xander, Thanks for both of you for the help provided throughout. I figured that out and its working now. I have tested few rasters manually.  ... View more

Dan, Please see my code complete code below which extract variable from netCDf file, project and save as .tif file.  import arcpy, sys from arcpy import env from arcpy.sa import * import time # Input data source arcpy.env.workspace = r"C:/Users/Desktop/weekly_avg_project" arcpy.env.overwriteOutput = True # Set output folder OutputFolder = r"C:/Users/Desktop/weekly_avg_project/output" # Loop through a list of files in the workspace NCfiles = arcpy.ListFiles("*.nc") start = time.time()‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍ NCfiles = arcpy.ListFiles("*.nc") start = time.time() for filename in NCfiles: print("Processing: " + filename) inNCfiles = arcpy.env.workspace + "//" + filename fileroot = filename[0:(len(filename)-3)] + '_var86' TempLayerFile = "SM_amount" outRaster = OutputFolder + "//" + fileroot # Process: Make NetCDF Raster Layer outfile = arcpy.MakeNetCDFRasterLayer_md(filename, "var86", "lon", "lat", TempLayerFile, "", "", "BY_VALUE") proj_ras = arcpy.ProjectRaster_management(outfile, "C:/Users/Desktop/weekly_avg_project/project_ras","PROJCS['USA_Contiguous_Albers_Equal_Area_Conic',GEOGCS['GCS_North_American_1983',DATUM['D_North_American_1983',SPHEROID['GRS_1980',6378137.0,298.257222101]],PRIMEM['Greenwich',0.0],UNIT['Degree',0.0174532925199433]],PROJECTION['Albers'],PARAMETER['False_Easting',0.0],PARAMETER['False_Northing',0.0],PARAMETER['Central_Meridian',-96.0],PARAMETER['Standard_Parallel_1',29.5],PARAMETER['Standard_Parallel_2',45.5],PARAMETER['Latitude_Of_Origin',37.5],UNIT['Meter',1.0]]", "BILINEAR", "10 10", "NAD_1983_To_WGS_1984_1", "", "GEOGCS['GCS_WGS_1984',DATUM['D_WGS_1984',SPHEROID['WGS_1984',6378137.0,298.257223563]],PRIMEM['Greenwich',0.0],UNIT['Degree',0.0174532925199433]]") # Process: Copy Raster arcpy.CopyRaster_management(proj_ras, outRaster + ".tif", "", "", "", "NONE", "NONE", "")‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍ ... View more

I used the above Geographic projection to be compatible with other shapefiles and rasters. But, in the code above I convert the projection from "GCS_North_American_1983" to "USA_Contiguous_Albers_Equal_Area_Conic projection" so it gives the cell size in meters.  Are you saying, extract and convert the netCDF variable to UTM projection? ... View more

No. The files did not have a projection defined. It was a variable extracted from netCDF files so I defined the projection and saved as a .tif raster. Please see the line that defines the raster below. proj_ras = arcpy.ProjectRaster_management(outfile, "C:/Users/Desktop/weekly_avg_project/project_ras","GEOGCS['GCS_North_American_1983',DATUM['D_North_American_1983',SPHEROID['GRS_1980',6378137.0,298.257222101]],PRIMEM['Greenwich',0.0],UNIT['Degree',0.0174532925199433]]", "NEAREST", "0.125 0.125", "NAD_1983_To_WGS_1984_1", "", "GEOGCS['GCS_WGS_1984',DATUM['D_WGS_1984',SPHEROID['WGS_1984',6378137.0,298.257223563]],PRIMEM['Greenwich',0.0],UNIT['Degree',0.0174532925199433]]")‍ ... View more

Yes. All the input rasters have a defined coordinate system.  ... View more

Dan, I have formatted the text above to better read the code. I have already defined the projection for all the rasters.  ... View more

Dan, I would like to take the above comment back. I am getting the following error message. I have defined the spatial reference for all the raster files in the "Top_SM" folder as well as for the "TownWilson" shape file.  Please see the error message and script below. Error: Executing: zonalMeanRev Start Time: Wed Jan 03 16:04:38 2018 Running script zonalMeanRev... Failed script zonalMeanRev... Traceback (most recent call last): File "C:\Users\Desktop\Arkansas Flood\My_Phyton_Scripts\zonal_mean_rev.py", line 39, in <module> prj_raster = arcpy.ProjectRaster_management(clip_raster, "#", "PROJCS['USA_Contiguous_Albers_Equal_Area_Conic',GEOGCS['GCS_North_American_1983',DATUM['D_North_American_1983',SPHEROID['GRS_1980',6378137.0,298.257222101]],PRIMEM['Greenwich',0.0],UNIT['Degree',0.0174532925199433]],PROJECTION['Albers'],PARAMETER['False_Easting',0.0],PARAMETER['False_Northing',0.0],PARAMETER['Central_Meridian',-96.0],PARAMETER['Standard_Parallel_1',29.5],PARAMETER['Standard_Parallel_2',45.5],PARAMETER['Latitude_Of_Origin',37.5],UNIT['Meter',1.0]]", "BILINEAR", "10 10", "", "", "GEOGCS['GCS_North_American_1983',DATUM['D_North_American_1983',SPHEROID['GRS_1980',6378137.0,298.257222101]],PRIMEM['Greenwich',0.0],UNIT['Degree',0.0174532925199433]]") File "c:\program files (x86)\arcgis\desktop10.4\arcpy\arcpy\management.py", line 9025, in ProjectRaster raise e ExecuteError: ERROR 999999: Error executing function. General function failure [week_4_avg_200609_var86_Clip5] General function failure [fras_blk_week_4_avg_200609_var86_Clip5] No spatial reference exists. The table was not found. [week_4_avg_200609_var86_Clip5] Failed to execute (ProjectRaster). Failed to execute (zonalMeanRev).‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍ ############################### import arcpy, os, arcinfo from arcpy import env from arcpy.sa import * arcpy.env.overwriteOutput = True arcpy.CheckOutExtension("Spatial") env.workspace = "H:/outRaster/Rasters/Top_SM" watershedFeat = "H:/Arkansas Flood/KS_new.gdb/TownWilson" outDir = "H:/Table/" x = arcpy.ListRasters() for raster in arcpy.ListRasters(): raster_name = os.path.basename(raster).rstrip(os.path.splitext(raster)[1]) ## Clip Raster Dataset by known extent - Left Bottom Right Top clip_raster = arcpy.Clip_management(raster, "-96.150 37.382844 -95.365 37.736391", "#", "#", "#", "NONE") # Process: Project Raster prj_raster = arcpy.ProjectRaster_management(clip_raster, "#", "PROJCS['USA_Contiguous_Albers_Equal_Area_Conic',GEOGCS['GCS_North_American_1983',DATUM['D_North_American_1983',SPHEROID['GRS_1980',6378137.0,298.257222101]],PRIMEM['Greenwich',0.0],UNIT['Degree',0.0174532925199433]],PROJECTION['Albers'],PARAMETER['False_Easting',0.0],PARAMETER['False_Northing',0.0],PARAMETER['Central_Meridian',-96.0],PARAMETER['Standard_Parallel_1',29.5],PARAMETER['Standard_Parallel_2',45.5],PARAMETER['Latitude_Of_Origin',37.5],UNIT['Meter',1.0]]", "BILINEAR", "10 10", "", "", "GEOGCS['GCS_North_American_1983',DATUM['D_North_American_1983',SPHEROID['GRS_1980',6378137.0,298.257222101]],PRIMEM['Greenwich',0.0],UNIT['Degree',0.0174532925199433]]") outTable = outDir + raster_name + "_TBL.dbf" arcpy.gp.ZonalStatisticsAsTable(watershedFeat,"T_R", prj_raster, outTable,"NODATA","MEAN") arcpy.CheckInExtension("Spatial")‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍ ... View more

Dan and Xander, I made it work and wrote a script to process over 500 rasters. Its working fine. Thanks a lot for both of you. ... View more

Dan and Xander, I have tried clipping the raster using the "Extract by Mask" tool in GIS for the grids shapefile (TownWilson) shown above. I made sure the raster and the shapefile are in the same Geographic coordinates projection. I have also set the environment variables as shown below. But, the output does not cover the grids. (See the images below)   Output Raster does not cover the grids. Appreciate if you can guide what I am doing wrong here.  Thanks a lot again.  ... View more

Xander, I am going to try "BILINEAR" using Resample tool. Cell size of the original raster is in degrees (0.125 x 0.125). So, I found the degrees equivalent to 10m resolution is 0.0000898271 and the set environment variables as well. (Please see the screen shots below.) Please let me know whether is this something you tried to obtain the Zonal Stat table above. Thanks, ... View more

I will try this. Any suggestion of what re-sampling technique I should use. Nearest/Bilinear/Cubic? The raster data has soil moisture.  ... View more

Thanks Xander. Did you re-sample raster to 10m resolution and then perform zonal statistics tool? What GIS tool you used to re-sample? I appreciate if you can let me know the steps followed to obtain the results/table above? ... View more

Thanks Dan. Do you suggest me to re-sample the raster for 10m resolution (as Xander finds out) before perform the zonal statistics.? ... View more