#! /usr/bin/env python

import numpy as np
import matplotlib.pyplot as plt

# Here we load the topography data file
ttt = '>i2' # MSB_INTEGER 16-bit
a = np.fromfile(open("megt90n000cb.img", "r"), dtype=ttt)
dt = np.dtype((ttt, (720, 1440)))

# Here we calibrate the size of the topography data file
x = np.squeeze(np.frombuffer(a,dt))/1000.

# Here we define the longitude/latitude grid
lat=720
lon=1440
latitude=np.zeros(lat,dtype='f')
longitude=np.zeros(lon,dtype='f')
for i in range(0,lon,1):
   longitude[i]=-180.+i/4.
for j in range(0,lat,1):
   latitude[j]=90.-j/4.
# Here we make the histogram
plt.figure(1)
plt.hist(x.flatten(),bins=30)
plt.ylabel('Number of bins')
plt.xlabel('elevation (km)')

# Here we plot the map of the topography
plt.figure(2)
plt.xlabel('longitude')
plt.ylabel('latitude')
plt.contourf(longitude,latitude,x)
plt.colorbar()
plt.show()
