# Python script to create a file with mountain peaks and their heights
# L. Fita, LMD March 2017
#
# To be used as example for PyNCplot `draw_ptZvals'
#
import numpy as np
from netCDF4 import Dataset as NetCDFFile
import generic_tools as gen
import nc_var_tools as ncvar

# A series of mountain peaks, their heights and their coordinates 
#   [DEGlon, MINlon, SEClon, DEGlat, MINlat, SEClat] (from Wikipedia)
peaks = {'Everest': [8848., 86., 55., 31., 27., 59., 17.],                           \
  'Mont Blanc': [4808.73, 6., 51., 53., 45., 49., 57.],                              \
  'Aneto': [3404., -0., -39., -28., 42., 37., 56.],                                  \
  'Aconcagua': [6962., -70., -0, -40., -32., -39., -12.],                            \
  'Denali': [6190., -151., -0., -27., 63., 04., 10.],                                \
  'Kilimanjaro': [5895., 37., 21, 35., -3.,-4.,-0.],                                 \
  'Ulrulu': [863., 131., 2., 11., -25., -20, -42.],                                  \
  'Elbrus': [5642., 42., 26., 21., 43., 21., 18.],                                   \
  'Vinson': [4892., -85., -32., -2., -78., -31., -31.],                              \
  'Enduwa Kombuglu': [4509., 145., 2., 0., -5., -48., -0.],                          \
  'Orizaba': [5636., -97., -16., -12., 19., 01., 48.] }

# Length of the peaks name
Lstring = 250

ofile = 'MountainPeaks.nc'

#######    ########
## MAIN
    #######

# Sorting peaks' names (just to be fancy...)
peaksnames = list(peaks.keys())
peaksnames.sort()

onc = NetCDFFile(ofile, 'w')

Nmountains = len(peaks.keys())

# Dimensions
newdim = onc.createDimension('mountain', Nmountains)
newdim = onc.createDimension('Lstring', Lstring)

# Variables
newvarname = onc.createVariable('name', 'c', ('mountain', 'Lstring'))
ncvar.basicvardef(newvarname,'name','Name of the peak','-')
newattr = ncvar.set_attribute(newvarname,'coordinates','lon lat')
newvals = ncvar.writing_str_nc(newvarname, peaksnames, Lstring)

newvarlon = onc.createVariable('lon', 'f8', ('mountain'))
ncvar.basicvardef(newvarlon,'longitude','Longitude','degrees_East')
newattr = ncvar.set_attribute(newvarlon,'axis','X')
newattr = ncvar.set_attribute(newvarlon,'_CoordinateAxisType','Lon')

newvarlat = onc.createVariable('lat', 'f8', ('mountain'))
ncvar.basicvardef(newvarlat,'latitude','Latitude','degrees_Noth')
newattr = ncvar.set_attribute(newvarlat,'axis','Y')
newattr = ncvar.set_attribute(newvarlat,'_CoordinateAxisType','Lat')

newvarhgt = onc.createVariable('height', 'f4', ('mountain'))
ncvar.basicvardef(newvarhgt,'height','Height above sea level','m')

for keyv in peaks.keys():
    peakv = peaks[keyv]
    ipeak = peaksnames.index(keyv)

    newvarhgt[ipeak] = peakv[0]
    newvarlon[ipeak] = peakv[1] + peakv[2]/60. + peakv[3]/3600.
    newvarlat[ipeak] = peakv[4] + peakv[5]/60. + peakv[6]/3600.

    onc.sync()

# Global values
addglob = ncvar.add_global_PyNCplot(onc, 'create_TopoValues.py', 'main', '0.1')

onc.sync()
onc.close()

print "Successful written of '" + ofile + "' !!"