<jupyter_start><jupyter_text>The elevation of Earth  Load the Python modulesWe will use the `numpy` module to deal with data as arrays and we will use `matplotlib` to visualize data. These modules are imported as `np` and `plt` respectively so that we can access their functions with those codes as shorthand. Importing these different modules and remembering these shorthand codes will feel a bit cumbersome at first, but you will get used to it as the semester progresses.<jupyter_code>import numpy as np
import matplotlib.pyplot as plt<jupyter_output><empty_output><jupyter_text>Load topography data Inspect np.loadtxtFirst we'll load the files of global topography that we will be analyzing using the ```np.loadtxt``` function. Before we use the function, let us have a look at its docstring. By typing the name of a function with a `?` after it, we see a description of what the function does, what parameters need to be passed to it and what additional parameters can be provided. Another way to get this same information is to press `shift+tab` after the function name.<jupyter_code>#np.loadtxt?
#plt.contourf?<jupyter_output><empty_output><jupyter_text>Load the filesThere are three files: an array of latitudes (size 540 x 1; ```etopo20lats.txt'```), an array of longitudes (size 1081 x 1; ```etopo20lons.txt```), and a gridded array of topography with units of meters (size 540 x 1081; ```etopo20data.txt```). We will load the files and assign them to the variables ```lats```, ```lons```, and ```topo_grid```. Variables are names we make up that can be assigned to different Python data types — in this case numpy arrays.These files are within the data folder that is a subfolder to the folder that this notebook is in. To load them, we need to provide the ```path``` which is a string that tells Python where the file is relative to the current directory. For the latitude file, that path is ```'data/etopo20lats.txt'```<jupyter_code>lats = np.loadtxt('data/etopo20lats.txt') 
lons = np.loadtxt('data/etopo20lons.txt')
topo_grid = np.loadtxt('data/etopo20data.txt')<jupyter_output><empty_output><jupyter_text>Inspect the dataThese data are a global model of topography (points above sea level) and bathymetry (points below sea level) generated by the National Geophysical Data Center (NGDC), an office of the National Oceanic and Atmospheric Administration (NOAA). We are using the 20 arcsecond version of the data (1/3 of a degree). Even higher resolution data are available: https://data.nodc.noaa.gov/cgi-bin/iso?id=gov.noaa.ngdc.mgg.dem:316It helps to see the data, so let's ```print``` the latitude array and see that the numpy array is indeed an array of latitudes.<jupyter_code>#print(lats)
#lats.shape<jupyter_output><empty_output><jupyter_text>The topo_grid array is more complicated than the array of latitudes. The `lats` array is 1-dimensional `[90, 89, 88,..]` while the `topo_grid` array is two dimensional with a elevation for every latitude and longitude value. We can determine the shape of the topo_grid array using `.shape`.<jupyter_code>topo_grid.shape
x=topo_grid.shape
#x
#x[0]<jupyter_output><empty_output><jupyter_text>Visualizing the data distributionLet's use ```matplotlib``` to visual the distribution of elevations on Earth. The elevations are within the `topo_grid` array. First, we need to manipulate the data a bit. We want to have a 1-dimensional list of all the (540 x 1081 = 583740) (583740) element array values to plot the distribution of elevations. We can use the numpy function ```reshape``` to do this. Instead of a grid, the data will all be in a single row after it is reshaped. To do so, we are applying a function (in this case ```reshape```) to an object (in this case the topo_grid object) and saving the result as a new numpy array which we will call topo_array.<jupyter_code>topo_array = topo_grid.reshape(583740)
x[0]*x[1]
topo_array<jupyter_output><empty_output><jupyter_text>Plot histogram of elevation dataNow we can plot a simple histogram with the `matplotlib` function `hist` that shows the number of elevation points within a specified number of bins. Remember that we imported ```matplotlib.pyplot as plt``` so to access the ```hist``` function of matplotlib we can use ```plt.hist```.<jupyter_code>plt.hist(topo_array,bins=50,label='elevation data')
plt.xlabel('elevation (m)')
plt.ylabel('number of points')
plt.title('histogram of global elevation')
plt.show()<jupyter_output><empty_output><jupyter_text>**_Discussion question:_** *How would you describe the shape of this distribution of elevations?*<jupyter_code>from IPython.display import YouTubeVideo
# Why Earth Has Two Levels | Hypsometric Curve
# Video credit: MinuteEarth
# Start at 3s, End at 2:25
YouTubeVideo('KOv3FGVmRcA?t=3', width="1000", height="500")<jupyter_output><empty_output><jupyter_text>Calculate maximum and minimum of global elevationNumpy has helpful functions that allow us to determine properties of data. In this case, let's investigate the minimum and maximum of the topo array (the values are in meters).<jupyter_code>topo_max = np.max(topo_array)
topo_min = np.min(topo_array)
print(topo_max)
print(topo_min)<jupyter_output>6228.8125
-9026.625<jupyter_text>**_Discussion question:_** *The summit of Mount Everest is 8,848 meters. Why is the maximum in this data set different than this?*  Calculate mean of global elevation**_Code for you to write_**Calculate the mean of ```topo_array``` using the function ```np.mean()``` and assign it to the variable ```topo_mean```. The code to do so should look like ```topo_mean = np.mean(topo_array)```Calculate the standard deviation of ```topo_array``` using the function ```np.std()``` and assign it to the variable ```topo_sigma```.Determine the number of elevation points in the ```topo_array``` using the function ```len()``` and assign it to the variable ```topo_n```.<jupyter_code>### Calculate the mean, standard deviation and number of points here

### Use the print() function to see the topo_mean, topo_sigma and topo_n
topo_mean=np.mean(topo_array)
topo_sigma=np.std(topo_array)
topo_n=len(topo_array)

print(f'The mean elevation is {topo_mean:.1f} with a standard deviation of {topo_sigma:.1f}')
print(f'There are a {topo_n:d} points in the sample')<jupyter_output>The mean elevation is -1894.7 with a standard deviation of 2653.8
There are a 583740 points in the sample<jupyter_text>Generate synthetic elevations for a hypothetical Planet NormUsing ```topo_mean```,  ```topo_sigma```, and ```topo_n``` you calculated, let's generate synthetic elevations for a planet with a normal distribution of elevation. A normal distribution is also known as a Gaussian or a bell-curve distribution and can describe a distribution of data with a mean and a standard deviation. Let's call this made-up planet ```Planet Norm```.The function `np.random.normal` generates random points from a normal distribution and we can specify that it has the same mean as Earth's elveation (`topo_mean`), the same standard deviation of `topo_sigma` and that as many points are drawn from the distribution as we have in the data set (`topo_n`).<jupyter_code>Planet_Norm_topo = np.random.normal(topo_mean,topo_sigma,topo_n)
plt.hist(topo_array,bins=50,label='Earth elevation data',zorder=100)
plt.hist(Planet_Norm_topo,bins=50,label='Planet Norm elevation data')
plt.xlabel('elevation (m)')
plt.ylabel('number of points')
plt.title('histogram of global elevation')
plt.legend()
plt.show()<jupyter_output><empty_output><jupyter_text>Plot Hypsometric CurveAnother useful method of visualizing the distribution of our data is a cumulative distribution function (CDF), which for topography is called a hypsometric curve. The value on the hypsometric curve represents the percentage of Earth's surface with an elevation higher than the elevation at the corresponding value.Let us sort the data by decreasing numbers using the ```sorted``` function. We can then use the ```np.linspace``` function to create a list from 0 to 100 that has the same number of points as our number of topography points.<jupyter_code>sorted_topo = sorted(topo_array,reverse=True)
#sorted_topo = np.asarray(sorted_topo)  #not needed
percent_topo = np.linspace(0,100,len(sorted_topo))<jupyter_output><empty_output><jupyter_text>Let use make a plot of sorted topography on the y-axis and the percent topography list on the x-axis.<jupyter_code>plt.plot(percent_topo,sorted_topo,label='Earth elevation')
plt.hlines(0,-1,101,label='sea level',linestyles='--')
plt.xlabel('Percent')
plt.ylabel('Elevation (m)')
plt.title('Hypsometric Curve')
plt.xlim(-1,101)
plt.legend(loc='best')
plt.show()
#plt.hlines?<jupyter_output><empty_output><jupyter_text>Plot the hyposmetric curve for Planet Norm**_Code for you to write_**Follow the same procedure as we did for Earth's elevation for our synthetic Planet Norm elevation to get it sorted. Be sure to give it a unique name from Earth's sorted topo (something like ```sorted_norm_topo```). Plot the Planet Norm curve on a figure along with Earth's hypsometric curve.<jupyter_code>## Generate the sort Planet Norm topo here<jupyter_output><empty_output><jupyter_text>**_Code for you to write_**We can observe where the Earth elevation and the Planet Norm hypsometric curves cross zero to determine what percentage of the planet's elevation is above sea level. Put the plotting code from above in a new cell and use the `plt.xlim()` and `plt.ylim()` functions to zoom-in on the curves and visualize where the curves cross sea level.<jupyter_code>### Change the axis bounds to determine where the curves cross 0
plt.plot(percent_topo,sorted_topo,label='Earth elevation')
plt.plot(...)
plt.hlines(0,-1,101,label='sea level',linestyles='--')
plt.xlabel('Percent')
plt.ylabel('Elevation (m)')
plt.title('Hypsometric Curve')
plt.legend(loc='best')
plt.show()<jupyter_output><empty_output><jupyter_text>**_Discussion questions_***Is there more or less land above sea level on Earth than one would expect for a normal distribution about Earth's mean elevation? Why is this the case?* *How does this compare to the statement made in the video clip?*  Making a topographic mapHistograms and CDFs are two useful 1D visualizations of our data, but what about the distribution of topography across the globe? Again, we'll need to manipulate the data a bit. We need grids of latitude and longitude in order to map the topography. Like `topo_grid` they will be (size 540 x 1081). `lat_grid` will have rows of constant latitude, and `lon_grid` will have columns of constant longtiude. To do this we'll call the numpy functions `repeat`, `reshape`, and `tile`.***First lets look at this with simpler cases***<jupyter_code>#Let's create an array (vector) and then a repeated matrix in which values change in the row direction
x=([1, 2, 3, 4, 5])
y=np.repeat(x,3)
z=y.reshape(len(x),3)
z
#Now lets create a matrix in which values change in the column direction
a=([6,7,8])
b=np.repeat(a,len(x))
c=(b.reshape(len(a),len(x))).transpose()
c
#For the column changing direction a short can be used
d=np.tile(a,(len(x),1))
d
# repeat the array of latitudes 1081 times
x=topo_grid.shape
latitude =  np.repeat(lats,x[1]) 

# reshape that (583740,1) element array to (540,1081
lat_grid = latitude.reshape(x[0],x[1]) 

# repeat the array of longtitudes 540 times as rows
lon_grid = np.tile(lons,(x[0],1)) 
x<jupyter_output><empty_output><jupyter_text>Import cartopyWe are going to use the module `cartopy` to make maps using its coordinate reference tools so let use import it.<jupyter_code>import cartopy.crs as ccrs<jupyter_output><empty_output><jupyter_text>Use the `contourf` function to make the mapUse the grids of latitude and longtitude to plot a Robinson projection map of topography using the `cartopy` function `contourf`. This function is *slow* so do not worry if it takes a while to run. Note the 'levels' kwarg controls the number of color steps in contouring. Larger values are slower to plot than smaller values, but smaller values have limited resolution.Take a looks this explaination for the 'transform' kwarg in contourf, which can be confusing: https://scitools.org.uk/cartopy/docs/latest/tutorials/understanding_transform.html<jupyter_code>#plt.contourf?
fig = plt.figure(figsize=(10,8))
ax = plt.axes(projection=ccrs.Robinson())
ax.set_global()
plt.contourf(lon_grid, lat_grid, topo_grid, levels=40, transform=ccrs.PlateCarree())

ax.coastlines()
ax.gridlines()

plt.colorbar(ax=ax,orientation='horizontal',label='elevation (m)')

plt.show()<jupyter_output><empty_output><jupyter_text>Visualizing the shallowest and deepest parts of the oceanLet us use the `Orthographic` projection (A GLOBE!) to look at the large underwater mountain range of the Mid-Atlantic ridge. **_Code for you to write_**The code block below needs the `plt.contourf` function to be completed. Use the version of the inputs that we used above as a guide. In addition, add the `vmax` parameter and set it equal to zero so that the colorbar only corresponds to elevations below sea level.<jupyter_code>fig = plt.figure(figsize=(15,10))
ax = plt.axes(projection=ccrs.Orthographic(central_longitude=-40,
    central_latitude=20))
ax.set_global()
plt.contourf(..., ..., ..., levels=..., transform=ccrs.PlateCarree(), vmax=0)

ax.coastlines()
ax.gridlines()

plt.colorbar(ax=ax,orientation='horizontal',label='elevation (m)')

plt.show()<jupyter_output><empty_output><jupyter_text>Let us use the `Orthographic` projection to look at the Marianas Trench. A globe centered at a longitude of 145ºE and 20ºN will work well. **_Code for you to write_**The code block below needs the `plt.axes` and `plt.contourf` function to be completed. For the `plt.axes`, use the `ccrs.Orthographic` centered at a longitude of 145ºE and 20ºN.<jupyter_code>fig = plt.figure(figsize=(15,10))
ax = plt.axes(projection=ccrs.Orthographic(central_longitude=...,
    central_latitude=...))
ax.set_global()
plt.contourf(..., ..., ..., levels=40, transform=ccrs.PlateCarree(), vmax=0)

ax.coastlines()
ax.gridlines()

plt.colorbar(ax=ax,orientation='horizontal',label='elevation (m)')

plt.show()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.09793668739400792, 'fraction_numerical': 0.016322781232334653, 'mean_word_length': 5.907271839921913, 'pattern_counts': {'":': 0, '<': 79, '<?xml version=': 0, '>': 79, 'https://': 2, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '20', 'dir': 'jupyter-structured-clean-dedup', 'id': '355347', 'license': 'BSD-3-Clause', 'n_tokens_mistral': 4606, 'n_tokens_neox': 4252, 'n_words': 1721, 'path': 'week01_inclass/W01_inclass_elevation.ipynb', 'repo_name': 'ds-connectors/EPS-88-FA21', 'size': '14152'}