#/bin/bash

site=daint

case ${site} in

   daint)
      
      function plot_ts {
      
         infile=${1}
         outfile=${2}
      
         outfile=$(echo ${outfile} | sed -e "s/.png//")
      
         cdo graph,device=png ${infile} ${outfile}

         mv ${outfile}.1.png ${outfile}.png
      
      }
      
      
      function plot_2d {
      
         infile=${1}
         outfile=${2}
      
         var=$(cdo showvar ${infile} | sed -e "s/^ *//")
         
         case ${var} in
            totprec)
                      min=""
                      max=""
                      clevs=',list=0;.5;1;1.5;2;2.5;3.5;4;5;15'
                      colour=",colour_min=brown,colour_max=blue,colour_triad=acw"
                      ;;
            temp2|tsurf)
                      min=""
                      max=""
                      clevs=",interval=3"
                      colour=",colour_min=violet,colour_max=red,colour_triad=cw"
                      ;;
            wind10)
                      min=""
                      max=""
                      clevs=""
                      colour=",colour_min=violet,colour_max=red,colour_triad=cw"
                      ;;
            *)        
                      min=""
                      max=""
                      clevs=""
                      colour=""
                      ;;
         esac
      
         outfile=$(echo ${outfile} | sed -e "s/_${var}.png//")
      
         cdo grfill,device=png${clevs}${min}${max}${colour} ${infile} ${outfile} ${min} ${max} ${cint}
      }
      ;;


      *)      # site is not daint
     
      function plot_ts {

         infile=${1}
         outfile=${2}



         var=$(cdo showvar ${infile} | sed -e "s/^ *//")

         cat << __EOF__ | sed -e "s/^ *//" > plot.py
         
            import xarray as xr
            import numpy as np
            import matplotlib.pyplot as plt
             
            ds=xr.open_dataset('${infile}')
             
            #-- create figure and axes instances
            dpi = 100
            fig = plt.figure(figsize=(1100/dpi, 1100/dpi), dpi=dpi)
            ax  = fig.add_axes([0.1,0.1,0.8,0.9])
             
            ds.${var}.isel(lon=0,lat=0).plot()
             
            plt.savefig('${outfile}', bbox_inches='tight', dpi=dpi)
__EOF__

      python plot.py

      }


      function plot_2d {  
         
         infile=${1}
         outfile=${2}

         var=$(cdo showvar ${infile} | sed -e "s/^ *//")

         case ${var} in
            totprec)
                      clevs=',[0,.25,.5,.75,1,1.25,1.5,1.75,2,2.25,2.5,3,3.5,4,4.5,5,15]'
                      colour=",cmap=plt.cm.terrain_r,extend='max'"
                      ;;
            temp2|tsurf)
                      clevs=",np.arange(-40,30,5)"
                      colour=",cmap=plt.cm.jet,extend='both'"
                      ;;
            wind10)
                      clevs=",np.arange(0,15,1.5)"
                      colour=",cmap=plt.cm.jet,extend='both'"
                      ;;
            *)
                      clevs=""
                      colour=""
                      ;;
         esac



         cat << __EOF__ | sed -e "s/^ *//" > plot.py
            from   mpl_toolkits.basemap import Basemap, cm
            import matplotlib.pyplot as plt
            from   netCDF4 import Dataset as open_ncfile
            import numpy as np
            
            #-- open netcdf file
            nc = open_ncfile('${infile}')
            
            #-- read variable
            var = nc.variables['${var}'][:,:]
            lat = nc.variables['lat'][:]
            lon = nc.variables['lon'][:]
            
            #-- create figure and axes instances
            dpi = 100
            fig = plt.figure(figsize=(1100/dpi, 1100/dpi), dpi=dpi)
            ax  = fig.add_axes([0.1,0.1,0.8,0.9])
            
            #-- create map
            map = Basemap(projection='cyl',llcrnrlat= -90.,urcrnrlat= 90.,\
                          resolution='c',  llcrnrlon=-180.,urcrnrlon=180.)
            
            #-- draw coastlines, state and country boundaries, edge of map
            map.drawcoastlines()
            map.drawstates()
            map.drawcountries()
            
            #-- create and draw meridians and parallels grid lines
            map.drawparallels(np.arange( -90., 90.,30.),labels=[1,0,0,0],fontsize=10)
            map.drawmeridians(np.arange(-180.,180.,30.),labels=[0,0,0,1],fontsize=10)
            
            #-- convert latitude/longitude values to plot x/y values
            x, y = map(*np.meshgrid(lon,lat))
            
            #-- contour levels
            clevs = np.arange(210,320,5)
            
            #-- draw filled contours
            cnplot = map.contourf(x,y,var${clevs}${colour})
            
            #-- add colorbar
            cbar = map.colorbar(cnplot,location='bottom',pad="10%")      #-- pad: distance between map and colorbar
            cbar.set_label('${unit}')                                    #-- add colorbar title string
            
            #-- add plot title
            plt.title('${var}')
            
            #-- display on screen
            #plt.show()
            
            #-- maximize and save file
            plt.savefig('${outfile}', bbox_inches='tight', dpi=dpi)
            
__EOF__

      # run python script

      python plot.py

      }

      ;;

esac
