<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Other modules for geoscientists
Nikolay Koldunov
koldunovn@gmail.com
This is part of Python for Geosciences notes.
=============
Some of the things will not work on ZMAW computers (Iris, Cartopy).
Iris
Iris seeks to provide a powerful, easy to use, and community-driven Python library for analysing and visualising meteorological and oceanographic data sets. Kind of Ferret replacement. Developed in the Met Office by group of 7 full time developers. There are more than 300 active python users in Met Office.
With Iris you can
Step1:  This is how iris cube look like
Step2:  We can perform different operations on cubes. For example create zonal mean
Step3:  Here we plot timesiries from one point
Step4:  Section along longitude
Step5:  Cartopy
Cartopy is a library providing cartographic tools for python.
Some of the key features of cartopy are
Step6:  We only change projection
Step7:  One of the things it was originally created for is to handle datelines properly. Below is an example form the post by Filipe Fernandes, that demonstrate this feature. He show how to plot Challenger Expedition track in Basemap and Cartopy.
First download the data
Step8:  Basemap version
Step9:  Cartopy version
Step10:  Pydap
Pydap is a pure Python library implementing the Data Access Protocol, also known as DODS or OPeNDAP. You can use Pydap as a client to access hundreds of scientific datasets in a transparent and efficient way through the internet; or as a server to easily distribute your data from a variety of formats.
Step11:  We are going to access sea ice data from CliSAP-Integrated Climate Data Center (ICDC)
Step12:  F2PY
The F2PY project is created to unify the efforts of supporting easy connection between Fortran and Python languages. Example below is from Using Python and FORTRAN with F2py.
Create FORTRAN file (use %%file instead of %%writefile if you on IPython < 1.0)
Step13:  Compile it with f2py
Step14:  Import resulting fib1.so as python library
Step15:  Read some auto generated documentation
Step16:  Use fib function
Step17:  netCDF4-python
netCDF4-python is advanced Python interface to the netCDF version 4 library.
Among other features it can read data from a multi-file netCDF dataset.
Let's download one more file from NCEP reanalysis data
Step18:  The air variable now have 2 years of data
Step19:  It also have very nice functions for dates processing
Step20:  We can convert our time values to the datetime format, that python can work with
Step21:  Note that we don't have to apply scale factor and add offset, it's done automatically.
<ASSISTANT_TASK:>
Python Code:
import iris
import iris.quickplot as qplt
temperature = iris.load_cube('air.sig995.2012.nc')
qplt.contourf(temperature[0,:,:])
gca().coastlines()
Explanation: Other modules for geoscientists
Nikolay Koldunov
koldunovn@gmail.com
This is part of Python for Geosciences notes.
=============
Some of the things will not work on ZMAW computers (Iris, Cartopy).
Iris
Iris seeks to provide a powerful, easy to use, and community-driven Python library for analysing and visualising meteorological and oceanographic data sets. Kind of Ferret replacement. Developed in the Met Office by group of 7 full time developers. There are more than 300 active python users in Met Office.
With Iris you can:
Use a single API to work on your data, irrespective of its original format.
Read and write (CF-)netCDF, GRIB, and PP files.
Easily produce graphs and maps via integration with matplotlib and cartopy.
Here we load data from netCDF file in to cube object and plot first time step. Note automatic unpacking, and the title.
End of explanation
print temperature
Explanation: This is how iris cube look like:
End of explanation
zonal_mean = temperature.collapsed('latitude', iris.analysis.MEAN)
qplt.contourf(zonal_mean)
Explanation: We can perform different operations on cubes. For example create zonal mean:
End of explanation
#Code is a bit more complicated in order to fix issue with dates formating
fig = figure()
qplt.plot(temperature[:,10,10])
fig.autofmt_xdate()
Explanation: Here we plot timesiries from one point:
End of explanation
qplt.plot(temperature[0,:,10])
Explanation: Section along longitude:
End of explanation
import cartopy.crs as ccrs
ax = plt.axes(projection=ccrs.PlateCarree())
qplt.contourf(temperature[0,:,:])
gca().coastlines()
Explanation: Cartopy
Cartopy is a library providing cartographic tools for python.
Some of the key features of cartopy are:
object oriented projection definitions
point, line, polygon and image transformations between projections
integration to expose advanced mapping in matplotlib with a simple and intuitive interface
Simple plot:
End of explanation
ax = plt.axes(projection=ccrs.Mollweide())
qplt.contourf(temperature[0,:,:])
gca().coastlines()
ax = plt.axes(projection=ccrs.Robinson())
qplt.contourf(temperature[0,:,:])
gca().coastlines()
fig = figure(figsize=(7,7))
ax = plt.axes(projection=ccrs.NorthPolarStereo())
ax.set_extent([0, 360, 50, 90], crs=ccrs.PlateCarree())
qplt.contourf(temperature[0,:,:])
gca().coastlines()
Explanation: We only change projection:
End of explanation
!wget https://raw.github.com/ocefpaf/ocefpaf.github.io/master/downloads/notebooks/data/challenger_path.csv
Explanation: One of the things it was originally created for is to handle datelines properly. Below is an example form the post by Filipe Fernandes, that demonstrate this feature. He show how to plot Challenger Expedition track in Basemap and Cartopy.
First download the data:
End of explanation
kw = dict(color='#FF9900', linestyle='-', linewidth=1.5)
lon, lat = np.loadtxt('./challenger_path.csv', delimiter=',', unpack=True)
from mpl_toolkits.basemap import Basemap
def make_basemap(projection='robin', figsize=(10, 5), resolution='c'):
    fig, ax = plt.subplots(figsize=figsize)
    m = Basemap(projection=projection, resolution=resolution,
                lon_0=0, ax=ax)
    m.drawcoastlines()
    m.fillcontinents(color='0.85')
    parallels = np.arange(-60, 90, 30.)
    meridians = np.arange(-360, 360, 60.)
    m.drawparallels(parallels, labels=[1, 0, 0, 0])
    m.drawmeridians(meridians, labels=[0, 0, 1, 0])
    return fig, m
fig, m = make_basemap()
_ = m.plot(*m(lon, lat), **kw)
Explanation: Basemap version:
End of explanation
import cartopy.crs as ccrs
import cartopy.feature as cfeature
def make_cartopy(projection=ccrs.Robinson(), figsize=(10, 5), resolution='110m'):
    fig, ax = plt.subplots(figsize=figsize, subplot_kw=dict(projection=projection))
    ax.set_global()
    ax.coastlines(resolution=resolution, color='k')
    gl = ax.gridlines(draw_labels=False)  # Only PlateCarree and Mercator plots are currently supported.
    ax.add_feature(cfeature.LAND, facecolor='0.75')
    return fig, ax
fig, ax = make_cartopy(projection=ccrs.Robinson(), resolution='110m')
_ = ax.plot(lon, lat, transform=ccrs.Geodetic(), **kw)
Explanation: Cartopy version:
End of explanation
from pydap.client import open_url
Explanation: Pydap
Pydap is a pure Python library implementing the Data Access Protocol, also known as DODS or OPeNDAP. You can use Pydap as a client to access hundreds of scientific datasets in a transparent and efficient way through the internet; or as a server to easily distribute your data from a variety of formats.
End of explanation
dataset = open_url("http://icdc.zmaw.de/thredds/dodsC/amsre_asi_nh_2011")
print dataset
ice = dataset['icecon']
ice.shape
ice.attributes
imshow(squeeze(ice[0,:,:])*0.10000000149011612)
colorbar()
Explanation: We are going to access sea ice data from CliSAP-Integrated Climate Data Center (ICDC)
End of explanation
%%writefile FIB1.F
C FILE: FIB1.F
      SUBROUTINE FIB(A,N)
C
C     CALCULATE FIRST N FIBONACCI NUMBERS
C
      INTEGER N
      REAL*8 A(N)
      DO I=1,N
         IF (I.EQ.1) THEN
            A(I) = 0.0D0
         ELSEIF (I.EQ.2) THEN
            A(I) = 1.0D0
         ELSE 
            A(I) = A(I-1) + A(I-2)
         ENDIF
      ENDDO
      END
C END FILE FIB1.F
Explanation: F2PY
The F2PY project is created to unify the efforts of supporting easy connection between Fortran and Python languages. Example below is from Using Python and FORTRAN with F2py.
Create FORTRAN file (use %%file instead of %%writefile if you on IPython < 1.0):
End of explanation
!f2py -c -m --fcompiler=gnu95 fib1 FIB1.F
Explanation: Compile it with f2py:
End of explanation
import fib1
Explanation: Import resulting fib1.so as python library:
End of explanation
print fib1.__doc__
print fib1.fib.__doc__
Explanation: Read some auto generated documentation:
End of explanation
import numpy as np
a=np.zeros(15,'d') 
fib1.fib(a) 
print a 
Explanation: Use fib function:
End of explanation
!wget ftp://ftp.cdc.noaa.gov/Datasets/ncep.reanalysis/surface/air.sig995.2011.nc
from netCDF4 import MFDataset
f = MFDataset('air.sig995.????.nc')
f.variables
air = f.variables['air']
time = f.variables['time']
lat = f.variables['lat']
Explanation: netCDF4-python
netCDF4-python is advanced Python interface to the netCDF version 4 library.
Among other features it can read data from a multi-file netCDF dataset.
Let's download one more file from NCEP reanalysis data:
End of explanation
air.shape
Explanation: The air variable now have 2 years of data:
End of explanation
from netCDF4 import num2date
time.units
Explanation: It also have very nice functions for dates processing:
End of explanation
time_conv = num2date(time, time.units)
time_conv
fig = figure()
plot(time_conv, air[:,10,10])
fig.autofmt_xdate()
Explanation: We can convert our time values to the datetime format, that python can work with:
End of explanation
contourf(lat[:],time_conv,air[:,:,10])
colorbar()
Explanation: Note that we don't have to apply scale factor and add offset, it's done automatically.
End of explanation
<END_TASK>