Wow, that's a pretty cool dataset.
The text file didn't have a proper header so it requires a bit more manipulation. If the header was just a single row with column names, it would be easier.
The code below converts the input text file into a feature class with only the surface nitrate, the deepest non-zero nitrate, and the depth of the deepest non-zero nitrate value.
Everything including the resulting shapefile is in the attached folder. Hopefully the comments in the code will help you understand what it does. A good way to play with it is to run it row by row (block by block) in something like PyScripter
import arcpy
import os
csv = r'C:\Nitrate_GridP_Annual\Nitrate_GridP_Annual.csv'
out_fc = r'C:\Nitrate_GridP_Annual\ni.shp'
# read the csv into memeory
rows = []
with open(csv, 'r') as f:
for row in f:
rows.append(row.strip())
# The first two rows are not very useful so let's drop them
rows.pop(0) # drop the first row
depths = rows.pop(0) # drop the second row, which is now the first one
depths = map(int, depths[depths.find(':')+1:].split(','))
# Get lat, lon, Nitrate near the surface, and last non-zero measurement
newrows = []
for row in rows:
rowdata = map(float, row.split(','))
lat, lon = rowdata[0], rowdata[1]
nitrate = rowdata[2:]
n_near_surface = nitrate.pop(0) # first is near the surface
# get the deepest non-zero value and its depth
n_deepest = 0.0
depth = 0
i = 0
for i in range(len(nitrate)):
n = nitrate
if n > 0:
n_deepest = n
d = depths[i+1]
newrow = [lat, lon, n_near_surface, n_deepest, d]
newrows.append(newrow)
# write output to feature class
out_path, out_name = os.path.split(out_fc)
sr = arcpy.SpatialReference(4326) # define coordinate system (WGS84)
fc = arcpy.management.CreateFeatureclass(out_path, out_name, "POINT", spatial_reference=sr).getOutput(0)
arcpy.management.AddField(fc, "N_0", "DOUBLE")
arcpy.management.AddField(fc, "N_DEEP", "DOUBLE")
arcpy.management.AddField(fc, "DEPTH", "DOUBLE")
with arcpy.da.InsertCursor(fc, ["SHAPE@", "N_0", "N_DEEP", "DEPTH"]) as ic:
for row in newrows:
pt = arcpy.Point(row[1], row[0]) # remember that lon is x, lat is y
newrow = [pt, row[2], row[3], row[4]]
ic.insertRow(newrow)
del row
del ic
print fc
Please check that the outputs are correct! I rushed a bit while writing this.
Hope this helps,
Filip.