Monday, 2 July 2018

Granularity of elevation in GMTED2010 data / understanding shapefiles


Following on from Records of shapefile do not contain every field, where @scruss provided some very helpful Python code.


I seem to have a fundamental inability to understand a shapefile (reader object).


Each Pyshp shape consists of a record plus a geometry (bounding box).


However, each record also has “north, south, west, east” items, which have a granularity of one degree.


And the min/max/mean elevations never change over that one degree. Every record in that one degree seems ot have the same elevation values, although each has a different bbox.



I need more detail, and had hoped/expected that those elevations would reflect the bbox.


Am I out of luck? I want to extract the highest elevation for each 0.01 degree square.


The field names are


: 
['ID',
'SOURCE_ORG',
'SOURCE',
'EL_SURFACE',
'NORTH',
'SOUTH',

'WEST',
'EAST',
'X_SRCE_RES',
'Y_SRCE_RES',
'HORZ_UNIT',
'COORD_SYS',
'HORZ_DATUM',
'VERT_DATUM',
'VERT_UNIT',
'MIN_ELEV',

'MAX_ELEV',
'MEAN_ELEV',
'SDEV_ELEV',
'PROD_DATE']

and here are the first two records, followed by their bboxes.


: 
[268,
'Univ of Bristol',
'Ant Radar and Laser Alt DEM',

'Reflective',
-55, <===== these four don't change for many records, but the booxes do
-90,
-180,
180,
1000.0,
1000.0,
'Meter',
'Polar Sterographic',
'WGS 84',

'WGS 84',
'Meter',
-82, <=== these are the elvations
4211, <=== which remain constant for the 1 dgeree square
2152.694,
1127.631,
'31Aug2011']


:

[269, <=== a new record number
'Univ of Bristol',
'Ant Radar and Laser Alt DEM',
'Reflective',
-55, <=== stil the same 1 degree square
-90,
-180,
180,
1000.0,
1000.0,

'Meter',
'Polar Sterographic',
'WGS 84',
'WGS 84',
'Meter',
-82, <=== still the same elevation data
4211,
2152.694,
1127.631,
'31Aug2011']


and the two bboxes :


[-103.71666666971798, -72.74206576893341, -103.58333333639, -72.71706576893436]

[-126.75833333546299, -73.17539910224943, -126.64166666880098, -73.14206576891738]

Here is the code that produces the above data


            data = shapefile.Reader(arguments.data_path + '\\' + fileName)

# get field names, skipping deleted flag

field_names = []
for f in data.fields[1:]:
field_names.append((f[0]))

shapeRecords = data.iterShapeRecords()

for shape in shapeRecords:
geometry = shape.shape
record = shape.record


No comments:

Post a Comment

arcpy - Changing output name when exporting data driven pages to JPG?

Is there a way to save the output JPG, changing the output file name to the page name, instead of page number? I mean changing the script fo...