| 3 | download required datset from net |
| 4 | |
| 5 | def Etopo(lon_area, lat_area, resolution): |
| 6 | Input |
| 7 | resolution: resolution of topography for both of longitude and latitude [deg] |
| 8 | (Original resolution is 0.0167 deg) |
| 9 | lon_area and lat_area: the region of the map which you want like [100, 130], [20, 25] |
| 10 | |
| 11 | Output |
| 12 | Mesh type longitude, latitude, and topography data |
| 13 | |
| 14 | |
| 15 | Read NetCDF data |
| 16 | data = Dataset("ETOPO1_Ice_g_gdal.grd", "r") |
| 17 | |
| 18 | |
| 19 | Get data |
| 20 | lon_range = data.variables['x_range'][:] |
| 21 | lat_range = data.variables['y_range'][:] |
| 22 | topo_range = data.variables['z_range'][:] |
| 23 | spacing = data.variables['spacing'][:] |
| 24 | dimension = data.variables['dimension'][:] |
| 25 | z = data.variables['z'][:] |
| 26 | lon_num = dimension[0] |
| 27 | lat_num = dimension[1] |
| 28 | |
| 29 | Prepare array |
| 30 | lon_input = np.zeros(lon_num); lat_input = np.zeros(lat_num) |
| 31 | for i in range(lon_num): |
| 32 | lon_input[i] = lon_range[0] + i * spacing[0] |
| 33 | for i in range(lat_num): |
| 34 | lat_input[i] = lat_range[0] + i * spacing[1] |
| 35 | |
| 36 | Create 2D array |
| 37 | lon, lat = np.meshgrid(lon_input, lat_input) |
| 38 | |
| 39 | Convert 2D array from 1D array for z value |
| 40 | topo = np.reshape(z, (lat_num, lon_num)) |
| 41 | |
| 42 | Skip the data for resolution |
| 43 | if ((resolution < spacing[0]) | (resolution < spacing[1])): |
| 44 | print('Set the highest resolution') |
| 45 | else: |
| 46 | skip = int(resolution/spacing[0]) |
| 47 | lon = lon[::skip,::skip] |