diff --git a/.gitignore b/.gitignore
index 78f64563e5784a02a66632b07f387220c7a9c16d..d4b7d99896683ab8a6d3531ac8c7e24e75a74d5f 100644
--- a/.gitignore
+++ b/.gitignore
@@ -5,7 +5,8 @@ __pycache__/
# C extensions
*.so
-
+*.nc
+*.xml
# Distribution / packaging
.Python
build/
@@ -14,6 +15,7 @@ dist/
downloads/
eggs/
.eggs/
+.idea/
lib/
lib64/
parts/
@@ -24,6 +26,8 @@ wheels/
.installed.cfg
*.egg
MANIFEST
+outputs
+.DS_Store
# PyInstaller
# Usually these files are written by a python script from a template
diff --git a/DESCRIPTION.rst b/DESCRIPTION.rst
new file mode 100644
index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391
diff --git a/README.md b/README.md
deleted file mode 100644
index 70c514e70958cf52c320848b2015b8ac4f3ce93f..0000000000000000000000000000000000000000
--- a/README.md
+++ /dev/null
@@ -1,20 +0,0 @@
-# PyNEMO
-
-To be udated soon. This work springboards from the [PyNEMO](http://pynemo.readthedocs.io/en/latest/index.html) Project.
-
-## What is this repository for? ##
-
-## How do I get set up? ##
-
-## Contribution guidelines ##
-
-## Who do I talk to? ##
-
-* Repo owner or admin
-
-jdha
-
-* Other community or team contact
-
-
-For more information regarding the use and development of PyNEMO see: [PyNEMO Wiki](https://github.com/jdha/PyNEMO/wiki)
diff --git a/README.rst b/README.rst
new file mode 100644
index 0000000000000000000000000000000000000000..af36ede7874795d9c6086ca75026080e38ed6b89
--- /dev/null
+++ b/README.rst
@@ -0,0 +1,114 @@
+PyNEMO
+======
+
+To be updated soon. This work springboards from the `PyNEMO Project <http://pynemo.readthedocs.io/en/latest/index.html/>`_.
+
+What is this repository for?
+----------------------------
+
+How do I get set up?
+--------------------
+
+Steps to take to install PyNEMO, creating a specific conda virtual environment is highly recommended.
+`click here for more about virtual enviroments <https://docs.conda.io/projects/conda/en/latest/user-guide/tasks/manage-environments.html/>`_
+
+- Install Git (outside scope of this readme)
+- Clone PyNEMO repository::
+
+ $ git clone https://github.com/NOC-MSM/PyNEMO.git
+
+- Install Conda, either Anaconda or Miniconda (outside scope of this readme)
+- Create conda environment for PyNEMO::
+
+ $ cd to/PyNEMO/directory
+ $ conda env create -f environment_pynemo.yml
+
+- Activate the new virtual environment::
+
+ $ source activate pynemo_env
+
+- Install Jave JRE (outside scope of this readme) and link libjvm.dylib to LD_LIBRARY_PATH variable::
+
+ $ export LD_LIBRARY_PATH=/path/to/java/library/folder/containing/libjvm.dylib:$LD_LIBARY_PATH # see notes below
+
+- Install PyNEMO::
+
+ $ cd /location/of/pynemo/repo
+ $ python setup.py build
+ $ python setup.py install
+
+This should result in PyNEMO being installed in the virtual environment, and can be checked by entering::
+
+ $ pynemo -v
+
+Resulting in a help usage prompt::
+
+ $ usage: pynemo -g -s <namelist.bdy>
+
+The virtual environment can be deactivated to return you to the normal prompt by typing::
+
+$ conda deactivate
+To reactivate, the following needs to be typed::
+
+ $ source activate pynemo_env
+
+
+To use PyNEMO, the following command is entered: (the example will run an benchmarking test)::
+
+ $ pynemo -s /path/to/namelist/file (e.g. PyNEMO/inputs/namelist_remote.bdy)
+
+**Additional NOTES**
+
+For Macbook Pro 2015, macOS Mojave and Java SDK 13 and JRE 8 the following path for the libjvm library should be correct::
+
+ /Library/Java/JavaVirtualMachines/jdk-13.0.1.jdk/Contents/Home/lib/server
+
+Resulting in the following command: (this will be different for different java versions and operating systems)::
+
+ $ export LD_LIBRARY_PATH=/Library/Java/JavaVirtualMachines/jdk-13.0.1.jdk/Contents/Home/lib/server:$LD_LIBRARY_PATH
+
+For an iMac 2013, macOS Catalina and JRE 8 only the followinng path was found to be correct::
+
+ /Library/Internet\ Plug-Ins/JavaAppletPlugin.plugin/Contents/Home/lib/server
+
+With the following command being required to set the environment variable::
+
+ $ export LD_LIBRARY_PATH=/Library/Internet\ Plug-Ins/JavaAppletPlugin.plugin/Contents/Home/lib/server:$LD_LIBRARY_PATH
+
+The conda environment creation command has not yet been tested. The yml document (can be opened using text editor) gives a list of all the modules and their versions that are required for PyNEMO so a environment can be constructed using this document as reference (or if you use pip!)
+
+**Update** conda environment yaml file has been tested (and works!) on a Macbook Pro 2015 and iMac 2013 running Anaconda 3.7 and Miniconda 3.7 respectively.
+
+Contribution guidelines
+-----------------------
+
+Bench Marking Tests
+-------------------
+
+The PyNEMO module can be tested using the bench marking namelist bdy file in the inputs folder. To check the outputs of the benchmark test, these can be visualised using the plotting script within the test_scripts folder. The following steps are required,
+
+- Run PyNEMO using the namelist file in the inputs folder (namelist_remote.bdy) e.g.::
+
+ $ pynemo -s /path/to/namelist/file
+
+- This will create two output files coordinates.bdy.nc and NNA_R12_bdyT_y1979)m11.nc in an outputs folder
+
+- To check the coordinates.bdy.nc has the correct boundary points, the script bdy_coords_plot.py will plot the domain boundaries and shown the different locations of the rim width (increasing number should go inwards) This script is located in the test_scripts folder.
+
+- The result should look like this (if using the current benchmark data)
+
+.. image:: /screenshots/example_bdy_coords.png
+ :width: 800
+ :alt: Example BDY coords output
+
+Who do I talk to?
+-----------------
+
+* Repo owner or admin
+
+ jdha
+
+* Other community or team contact
+
+
+For more information regarding the use and development of PyNEMO see: [PyNEMO Wiki](https://github.com/jdha/PyNEMO/wiki)
diff --git a/README_markdown.md b/README_markdown.md
new file mode 100644
index 0000000000000000000000000000000000000000..3ad308b12af6f4a9f4b5753e72f660cc1697281c
--- /dev/null
+++ b/README_markdown.md
@@ -0,0 +1,104 @@
+# PyNEMO
+
+To be udated soon. This work springboards from the [PyNEMO](http://pynemo.readthedocs.io/en/latest/index.html) Project.
+
+## What is this repository for? ##
+
+## How do I get set up? ##
+
+Steps to take to install PyNEMO, creating a specific conda virtual environment is highly recommended. [click here for more about virtual enviroments](https://docs.conda.io/projects/conda/en/latest/user-guide/tasks/manage-environments.html)
+
+- Install Conda, either Anaconda or Miniconda (outside of this readme)
+- Create conda environment for PyNEMO
+```
+$ cd to/PyNEMO/directory
+```
+```
+$ conda env create -f environment_pynemo.yml
+```
+- Activate the new virtual environment
+```
+$ source activate pynemo_env
+```
+- Install Jave JRE (outside this readme) and link libjvm.dylib to LD_LIBRARY_PATH variable
+```
+$ export LD_LIBRARY_PATH=/path/to/java/library/folder/containing/libjvm.dylib:$LD_LIBARY_PATH # see notes below
+```
+- Install Git (outside this readme)
+```
+$ git clone https://github.com/NOC-MSM/PyNEMO.git
+```
+- Install PyNEMO:
+```
+$ cd /location/of/pynemo/repo
+```
+```
+$ python setup.py build
+```
+```
+$ python setup.py install
+```
+
+This should result in PyNEMO being installed in the virtual environment, and can be checked by entering:
+```
+$ pynemo -v
+```
+
+Resulting in a help usage prompt:
+```
+$ usage: pynemo -g -s <namelist.bdy>
+```
+
+The virtual environment can be deactivated to return you to the normal prompt by typing:
+```
+$ conda deactivate
+```
+
+To reactivate, the following needs to be typed
+```
+$ source activate pynemo_env
+```
+
+To use PyNEMO, the following command is entered: (the example will run an benchmarking test)
+```
+$ pynemo -s /path/to/namelist/file (e.g. PyNEMO/inputs/namelist_remote.bdy)
+```
+
+**Additional NOTES**
+
+for MacOs and Java SDK 13 and JRE 8 the following path should be correct - /Library/Java/JavaVirtualMachines/jdk-13.0.1.jdk/Contents/Home/lib/server
+
+Resulting in the following command: (this will be different for different java versions and operating systems)
+```
+$ export LD_LIBRARY_PATH=/Library/Java/JavaVirtualMachines/jdk-13.0.1.jdk/Contents/Home/lib/server:$LD_LIBRARY_PATH
+```
+The conda environment creation command has not yet been tested. The yml document (can be opened using text editor) gives a list of all the modules and their versions that are required for PyNEMO so a environment can be constructed using this document as reference (or if you use pip!)
+
+## Contribution guidelines ##
+
+## Bench Marking Tests ##
+
+The PyNEMO module can be tested using the bench marking namelist bdy file in the inputs folder. To check the outputs of the benchmark test, these can be visualised using the plotting script within the test_scripts folder. The following steps are required,
+
+- Run PyNEMO using the namelist file in the inputs folder (namelist_remote.bdy) e.g.
+
+ - $ pynemo -s /path/to/namelist/file
+
+- This will create two output files coordinates.bdy.nc and NNA_R12_bdyT_y1979)m11.nc in an outputs folder
+
+- To check the coordinates.bdy.nc has the correct boundary points, the script bdy_coords_plot.py will plot the domain boundaries and shown the different locations of the rim width (increasing number should go inwards) This script is located in the test_scripts folder.
+
+- The result should look like this (if using the current benchmark data)
+
+
+
+## Who do I talk to? ##
+
+* Repo owner or admin
+
+jdha
+
+* Other community or team contact
+
+
+For more information regarding the use and development of PyNEMO see: [PyNEMO Wiki](https://github.com/jdha/PyNEMO/wiki)
diff --git a/docs/.nojekyll b/docs/.nojekyll
new file mode 100644
index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391
diff --git a/environment_pynemo.yml b/environment_pynemo.yml
new file mode 100644
index 0000000000000000000000000000000000000000..e8df8c86ca996ad9d9c2003724730cdaed585840
--- /dev/null
+++ b/environment_pynemo.yml
@@ -0,0 +1,181 @@
+name: pynemo_env
+channels:
+ - conda-forge
+ - anaconda
+ - srikanthnagella
+ - defaults
+dependencies:
+ - alabaster=0.7.12=py27_0
+ - appnope=0.1.0=py27hb466136_0
+ - asn1crypto=1.0.1=py27_0
+ - astroid=1.6.5=py27_0
+ - attrs=19.2.0=py_0
+ - babel=2.7.0=py_0
+ - backports=1.0=py_2
+ - backports.functools_lru_cache=1.5=py_2
+ - backports.shutil_get_terminal_size=1.0.0=py27_2
+ - backports_abc=0.5=py27h6972548_0
+ - basemap=1.2.0=py27h0acbc05_0
+ - blas=1.0=mkl
+ - bleach=3.1.0=py27_0
+ - bzip2=1.0.8=h1de35cc_0
+ - ca-certificates=2019.10.16=0
+ - cartopy=0.16.0=py27h9263bd1_0
+ - certifi=2019.9.11=py27_0
+ - cffi=1.12.3=py27hb5b8e2f_0
+ - cftime=1.0.3.4=py27h1d22016_1001
+ - chardet=3.0.4=py27_1003
+ - cloudpickle=1.2.2=py_0
+ - configparser=4.0.2=py27_0
+ - cryptography=2.7=py27ha12b0ac_0
+ - curl=7.65.3=ha441bb4_0
+ - cycler=0.10.0=py27hfc73c78_0
+ - cython=0.29.13=py27h0a44026_0
+ - dbus=1.13.6=h90a0687_0
+ - decorator=4.4.0=py27_1
+ - defusedxml=0.6.0=py_0
+ - docutils=0.15.2=py27_0
+ - entrypoints=0.3=py27_0
+ - enum34=1.1.6=py27_1
+ - expat=2.2.6=h0a44026_0
+ - freetype=2.9.1=hb4e5f40_0
+ - functools32=3.2.3.2=py27_1
+ - futures=3.3.0=py27_0
+ - geos=3.6.2=h5470d99_2
+ - gettext=0.19.8.1=h15daf44_3
+ - glib=2.56.2=hd9629dc_0
+ - hdf4=4.2.13=h39711bb_2
+ - hdf5=1.10.1=ha036c08_1
+ - icu=58.2=h4b95b61_1
+ - idna=2.8=py27_0
+ - imagesize=1.1.0=py27_0
+ - intel-openmp=2019.4=233
+ - ipaddress=1.0.22=py27_0
+ - ipykernel=4.10.0=py27_0
+ - ipython=5.8.0=py27_0
+ - ipython_genutils=0.2.0=py27h8b9a179_0
+ - ipywidgets=7.5.1=py_0
+ - isort=4.3.21=py27_0
+ - jedi=0.15.1=py27_0
+ - jinja2=2.10.3=py_0
+ - jpeg=9b=he5867d9_2
+ - jsonschema=3.0.2=py27_0
+ - jupyter=1.0.0=py27_7
+ - jupyter_client=5.3.3=py27_1
+ - jupyter_console=5.2.0=py27_1
+ - jupyter_core=4.5.0=py_0
+ - keyring=18.0.0=py27_0
+ - kiwisolver=1.1.0=py27h0a44026_0
+ - krb5=1.16.1=hddcf347_7
+ - lazy-object-proxy=1.4.2=py27h1de35cc_0
+ - libcurl=7.65.3=h051b688_0
+ - libcxx=4.0.1=hcfea43d_1
+ - libcxxabi=4.0.1=hcfea43d_1
+ - libedit=3.1.20181209=hb402a30_0
+ - libffi=3.2.1=h475c297_4
+ - libgfortran=3.0.1=h93005f0_2
+ - libiconv=1.15=hdd342a3_7
+ - libnetcdf=4.5.0=h42fd751_7
+ - libpng=1.6.37=ha441bb4_0
+ - libsodium=1.0.16=h3efe00b_0
+ - libssh2=1.8.2=ha12b0ac_0
+ - libtiff=4.0.10=hcb84e12_2
+ - libxml2=2.9.9=hf6e021a_1
+ - libxslt=1.1.33=h33a18ac_0
+ - lxml=4.4.1=py27hef8c89e_0
+ - markupsafe=1.1.1=py27h1de35cc_0
+ - matplotlib=2.2.3=py27h54f8f79_0
+ - mccabe=0.6.1=py27_1
+ - mistune=0.8.4=py27h1de35cc_0
+ - mkl=2019.4=233
+ - mkl-service=2.3.0=py27hfbe908c_0
+ - mkl_fft=1.0.14=py27h5e564d8_0
+ - mkl_random=1.1.0=py27ha771720_0
+ - nbconvert=5.6.0=py27_1
+ - nbformat=4.4.0=py27hddc86d0_0
+ - ncurses=6.1=h0a44026_1
+ - netcdf4=1.3.1=py27he3ffdca_2
+ - notebook=5.7.8=py27_0
+ - numpy=1.16.5=py27hacdab7b_0
+ - numpy-base=1.16.5=py27h6575580_0
+ - numpydoc=0.9.1=py_0
+ - olefile=0.46=py27_0
+ - openssl=1.1.1d=h1de35cc_3
+ - owslib=0.18.0=py_0
+ - packaging=19.2=py_0
+ - pandoc=2.2.3.2=0
+ - pandocfilters=1.4.2=py27_1
+ - parso=0.5.1=py_0
+ - pathlib2=2.3.5=py27_0
+ - pcre=8.43=h0a44026_0
+ - pexpect=4.7.0=py27_0
+ - pickleshare=0.7.5=py27_0
+ - pillow=6.2.0=py27hb68e598_0
+ - pip=19.2.3=py27_0
+ - proj4=5.0.1=h1de35cc_0
+ - prometheus_client=0.7.1=py_0
+ - prompt_toolkit=1.0.15=py27h4a7b9c2_0
+ - psutil=5.6.3=py27h1de35cc_0
+ - ptyprocess=0.6.0=py27_0
+ - pycodestyle=2.5.0=py27_0
+ - pycparser=2.19=py27_0
+ - pyepsg=0.4.0=py27_0
+ - pyflakes=2.1.1=py27_0
+ - pygments=2.4.2=py_0
+ - pyjnius=1.4=py27_0
+ - pylint=1.9.2=py27_0
+ - pyopenssl=19.0.0=py27_0
+ - pyparsing=2.4.2=py_0
+ - pyproj=1.9.5.1=py27h833a5d7_1
+ - pyqt=4.11.4=py27_4
+ - pyrsistent=0.15.4=py27h1de35cc_0
+ - pyshp=2.1.0=py_0
+ - pysocks=1.7.1=py27_0
+ - python=2.7.16=h97142e2_7
+ - python-dateutil=2.8.0=py27_0
+ - python.app=2=py27_9
+ - pytz=2019.3=py_0
+ - pyzmq=18.1.0=py27h0a44026_0
+ - qt=4.8.7=1
+ - qtawesome=0.6.0=py_0
+ - qtconsole=4.5.5=py_0
+ - qtpy=1.9.0=py_0
+ - readline=7.0=h1de35cc_5
+ - requests=2.22.0=py27_0
+ - rope=0.14.0=py_0
+ - scandir=1.10.0=py27h1de35cc_0
+ - scipy=1.2.1=py27h1410ff5_0
+ - seawater=3.3.4=py_1
+ - send2trash=1.5.0=py27_0
+ - setuptools=41.4.0=py27_0
+ - shapely=1.6.4=py27h20de77a_0
+ - simplegeneric=0.8.1=py27_2
+ - singledispatch=3.4.0.3=py27he22c18d_0
+ - sip=4.18=py27_0
+ - six=1.12.0=py27_0
+ - snowballstemmer=2.0.0=py_0
+ - sphinx=1.8.5=py27_0
+ - sphinxcontrib=1.0=py27_1
+ - sphinxcontrib-websupport=1.1.2=py_0
+ - spyder=3.2.8=py27_0
+ - spyder-kernels=1.4.0=py27_0
+ - sqlite=3.30.0=ha441bb4_0
+ - subprocess32=3.5.4=py27h1de35cc_0
+ - terminado=0.8.2=py27_0
+ - testpath=0.4.2=py27_0
+ - thredds_crawler=1.0.0=py27_0
+ - tk=8.6.8=ha441bb4_0
+ - tornado=5.1.1=py27h1de35cc_0
+ - traitlets=4.3.3=py27_0
+ - typing=3.7.4.1=py27_0
+ - urllib3=1.24.2=py27_0
+ - wcwidth=0.1.7=py27h817c265_0
+ - webencodings=0.5.1=py27_1
+ - wheel=0.33.6=py27_0
+ - widgetsnbextension=3.5.1=py27_0
+ - wrapt=1.11.2=py27h1de35cc_0
+ - wurlitzer=1.0.3=py27_0
+ - xz=5.2.4=h1de35cc_4
+ - zeromq=4.3.1=h0a44026_3
+ - zlib=1.2.11=h1de35cc_3
+ - zstd=1.3.7=h5bba6e5_0
diff --git a/inputs/namelist_remote.bdy b/inputs/namelist_remote.bdy
index 09d77854d09e8d293fde27ab468988ba779fa7e7..72b58c848840f7985bcd296a2fc4f6b475442701 100644
--- a/inputs/namelist_remote.bdy
+++ b/inputs/namelist_remote.bdy
@@ -30,18 +30,18 @@
!------------------------------------------------------------------------------
! grid information
!------------------------------------------------------------------------------
- sn_src_hgr = './benchmark/grid_low_res_C/mesh_hgr.nc'
- sn_src_zgr = './benchmark/grid_low_res_C/mesh_zgr.nc'
- sn_dst_hgr = './benchmark/grid_C/mesh_hgr_zps.nc'
- sn_dst_zgr = './benchmark/grid_C/mesh_zgr_zps.nc'
- sn_src_msk = './benchmark/grid_low_res_C/mask.nc'
- sn_bathy = './benchmark/grid_C/NNA_R12_bathy_meter_bench.nc'
+ sn_src_hgr = 'http://opendap4gws.jasmin.ac.uk/thredds/noc_msm/dodsC/pynemo_grid_low_res_C/mesh_hgr.nc'
+ sn_src_zgr = 'http://opendap4gws.jasmin.ac.uk/thredds/noc_msm/dodsC/pynemo_grid_low_res_C/mesh_zgr.nc'
+ sn_dst_hgr = 'http://opendap4gws.jasmin.ac.uk/thredds/noc_msm/dodsC/pynemo_grid_C/mesh_hgr_zps.nc'
+ sn_dst_zgr = 'http://opendap4gws.jasmin.ac.uk/thredds/noc_msm/dodsC/pynemo_grid_C/mesh_zgr_zps.nc'
+ sn_src_msk = 'http://opendap4gws.jasmin.ac.uk/thredds/noc_msm/dodsC/pynemo_grid_low_res_C/mask.nc'
+ sn_bathy = 'http://opendap4gws.jasmin.ac.uk/thredds/noc_msm/dodsC/pynemo_grid_C/NNA_R12_bathy_meter_bench.nc'
!------------------------------------------------------------------------------
! I/O
!------------------------------------------------------------------------------
- sn_src_dir = '/Users/jdha/Projects/GitHub/PyNEMO/inputs/src_data.ncml' ! src_files/'
- sn_dst_dir = './outputs'
+ sn_src_dir = '/Users/thopri/Projects/PyNEMO/inputs/src_data_remote.ncml' ! src_files/'
+ sn_dst_dir = '/Users/thopri/Projects/PyNEMO/outputs'
sn_fn = 'NNA_R12' ! prefix for output files
nn_fv = -1e20 ! set fill value for output files
nn_src_time_adj = 0 ! src time adjustment
diff --git a/pynemo/nemo_bdy_zgrv2.py b/pynemo/nemo_bdy_zgrv2.py
index e91fb8ca9d7120b85dafc6c006d0ce81cbc7100f..92bf9edcfb7cf960fb6dbbff59f94f9be98f6437 100644
--- a/pynemo/nemo_bdy_zgrv2.py
+++ b/pynemo/nemo_bdy_zgrv2.py
@@ -68,6 +68,7 @@ class Depth:
v_ind = sub2ind(mbathy.shape, bdy_v[:,0], bdy_v[:,1])
v_ind2 = sub2ind(mbathy.shape, bdy_v[:,0], bdy_v[:,1] + 1)
+ [tmp_zt, tmp_zw] = e3_to_depth(np.squeeze(nc['e3t'][:,:,:,:]), np.squeeze(nc['e3w'][:,:,:,:]), nz)
# This is very slow
self.logger.debug( 'starting nc reads loop' )
for k in range(nz):
@@ -85,7 +86,8 @@ class Depth:
# jelt: replace 'load gdep[wt] with load e3[tw] and compute gdep[tw]
#wrk1 = nc['gdept'][0,k,:,:]#nc.variables['gdept'][0,k,:,:]
#wrk2 = nc['gdepw'][0,k,:,:]#nc.variables['gdepw'][0,k,:,:]
- [wrk1, wrk2] = e3_to_depth(nc['e3t'][0,k,:,:], nc['e3w'][0,k,:,:], nz)
+ #print 'e3t shape: ', nc['e3t_0'][:].shape
+ [wrk1, wrk2] = tmp_zt[k,:,:], tmp_zw[k,:,:]
# Replace deep levels that are not used with NaN
wrk2[mbathy + 1 < k + 1] = np.NaN
diff --git a/pynemo/reader/ncml.py b/pynemo/reader/ncml.py
index 804f1f65f8a585428d8b758340dc4a5f705d9a3a..5d1e4504f79954c920a5c57b174a9d98e3a8d86e 100644
--- a/pynemo/reader/ncml.py
+++ b/pynemo/reader/ncml.py
@@ -225,7 +225,7 @@ class Variable(object):
retval = data_array.copyToNDJavaArray()
#TODO: copy it into numpy instead of Java array and then convert to numpy
# convert to numpy array
- retval = np.asarray(retval)
+ retval = np.asarray(retval, dtype='float')
self.logger.info(retval.shape)
if np_input: #if an array is passed as selection
ret_dim_list = ()
diff --git a/pynemo/tests/bdy_coords.py b/pynemo/tests/bdy_coords.py
new file mode 100755
index 0000000000000000000000000000000000000000..5f801fffa9bbd0b20850bc610c951db156828237
--- /dev/null
+++ b/pynemo/tests/bdy_coords.py
@@ -0,0 +1,157 @@
+'''
+###
+## This is a subset of pynemo for Robinson. It should be able to read a mask
+## file and produce a coordinates.bdy.nc file (jdha@noc.ac.uk)
+###
+'''
+
+# pylint: disable=E1103
+# pylint: disable=no-name-in-module
+
+#External imports
+from time import clock
+import numpy as np
+import logging
+
+#local imports
+from pynemo import nemo_bdy_setup as setup
+from pynemo import nemo_bdy_gen_c as gen_grid
+from pynemo import nemo_coord_gen_pop as coord
+
+from pynemo import nemo_bdy_source_coord as source_coord
+from pynemo import nemo_bdy_dst_coord as dst_coord
+from pynemo import nemo_bdy_ice
+from pynemo import pynemo_settings_editor
+from pynemo.utils import Constants
+
+from pynemo.gui.nemo_bdy_mask import Mask as Mask_File
+from PyQt4.QtGui import QMessageBox
+#import pickle
+
+
+logger = logging.getLogger(__name__)
+def process_bdy(setup_filepath=0, mask_gui=False):
+ """ Main entry to the processing of the bdy
+ Keyword arguments:
+ setup_filepath -- file path to bdy file
+ mask_gui -- whether gui to select the mask file needs to be poped up
+ """
+ #Logger
+ logger.info('START')
+ start = clock()
+ SourceCoord = source_coord.SourceCoord()
+ DstCoord = dst_coord.DstCoord()
+
+ logger.info(clock() - start)
+ start = clock()
+ Setup = setup.Setup(setup_filepath) # default settings file
+ settings = Setup.settings
+
+ logger.info(clock() - start)
+ ice = settings['ice']
+ logger.info('ice = %s', ice)
+
+ logger.info('Done Setup')
+ # default file, region settingas
+
+ start = clock()
+ bdy_msk = _get_mask(Setup, mask_gui)
+ logger.info(clock() - start)
+ logger.info('Done Mask')
+
+ DstCoord.bdy_msk = bdy_msk == 1
+ reload(gen_grid)
+ start = clock()
+ logger.info('start bdy_t')
+ grid_t = gen_grid.Boundary(bdy_msk, settings, 't')
+
+ logger.info(clock() - start)
+ start = clock()
+ logger.info('start bdy_u')
+ grid_u = gen_grid.Boundary(bdy_msk, settings, 'u')
+ logger.info('start bdy_v')
+
+ logger.info(clock() - start)
+ start = clock()
+ grid_v = gen_grid.Boundary(bdy_msk, settings, 'v')
+ logger.info('start bdy_f')
+
+ logger.info(clock() - start)
+ start = clock()
+
+ grid_f = gen_grid.Boundary(bdy_msk, settings, 'f')
+ logger.info('done bdy t,u,v,f')
+
+ logger.info(clock() - start)
+ start = clock()
+ if ice:
+ grid_ice = nemo_bdy_ice.BoundaryIce()
+ grid_ice.grid_type = 't'
+ grid_ice.bdy_r = grid_t.bdy_r
+
+ bdy_ind = {'t': grid_t, 'u': grid_u, 'v': grid_v, 'f': grid_f}
+
+ for k in bdy_ind.keys():
+ logger.info('bdy_ind %s %s %s', k, bdy_ind[k].bdy_i.shape, bdy_ind[k].bdy_r.shape)
+
+ start = clock()
+ co_set = coord.Coord(settings['dst_dir']+'/coordinates.bdy.nc', bdy_ind)
+ logger.info('done coord gen')
+ logger.info(clock() - start)
+ start = clock()
+ logger.info(settings['dst_hgr'])
+ co_set.populate(settings['dst_hgr'])
+ logger.info('done coord pop')
+
+ logger.info(clock() - start)
+
+ # may need to rethink grid info
+ # tracer 3d frs over rw
+ # tracer 2d frs over rw (i.e. ice)
+ # dyn 2d over 1st rim of T grid (i.e. ssh)
+ # dyn 2d over 1st rim
+ # dyn 2d frs over rw
+ # dyn 3d over 1st rim
+ # dyn 3d frs over rw
+
+def _get_mask(Setup, mask_gui):
+ """ This method reads the mask information from the netcdf file or opens a gui
+ to create a mask depending on the mask_gui input. return the mask data. The default mask
+ data is using bathymetry and applying a 1px halo
+ Keyword arguments:
+ Setup -- settings for bdy
+ mask_gui -- boolean to open mask gui.
+ """
+ bdy_msk = None
+ if mask_gui:
+ #Open the gui to create a mask
+ _, mask = pynemo_settings_editor.open_settings_dialog(Setup)
+ bdy_msk = mask.data
+ Setup.refresh()
+ else:
+ try:
+ #mask filename and mask file flag is set
+ if Setup.bool_settings['mask_file'] and Setup.settings['mask_file'] is not None:
+ mask = Mask_File(mask_file=Setup.settings['mask_file'])
+ bdy_msk = mask.data
+ elif Setup.bool_settings['mask_file']:
+ logger.error("Mask file is not given")
+ return
+ else: #no mask file specified then use default 1px halo mask
+ logger.warning("Using default mask with bathymetry!!!!")
+ mask = Mask_File(Setup.settings['bathy'])
+ mask.apply_border_mask(Constants.DEFAULT_MASK_PIXELS)
+ bdy_msk = mask.data
+ except ValueError: # why is this except here? as there is an else: statement TODO
+ print 'something wrong?'
+ return
+ if np.amin(bdy_msk) == 0:
+ # Mask is not set throw a warning message and set border to 1px.
+ logger.warning("Setting the mask to 1px border")
+ QMessageBox.warning(None,"pyNEMO", "Mask is not set, setting a 1 pixel border mask")
+ if bdy_msk is not None and 1 < bdy_msk.shape[0] and 1 < bdy_msk.shape[1]:
+ tmp = np.ones(bdy_msk.shape, dtype=bool)
+ tmp[1:-1, 1:-1] = False
+ bdy_msk[tmp] = -1
+
+ return bdy_msk
diff --git a/screenshots/example_bdy_coords.png b/screenshots/example_bdy_coords.png
new file mode 100644
index 0000000000000000000000000000000000000000..cd83e138a1bc0d264fabd8c3c0e25c2686696aa2
Binary files /dev/null and b/screenshots/example_bdy_coords.png differ
diff --git a/test_scripts/__init__.py b/test_scripts/__init__.py
new file mode 100644
index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391
diff --git a/test_scripts/bdy_coords_plot.py b/test_scripts/bdy_coords_plot.py
new file mode 100755
index 0000000000000000000000000000000000000000..9647cd5c41716e30eddb69b7774809abd7717e68
--- /dev/null
+++ b/test_scripts/bdy_coords_plot.py
@@ -0,0 +1,63 @@
+#!/usr/bin/env python2
+# -*- coding: utf-8 -*-
+"""
+Created on Thu Oct 24 11:28:04 2019
+
+@author: thopri
+"""
+
+# Hard file/folder paths in namelist file and in process_bdy call below need to be updated to suit user
+
+import matplotlib.pyplot as plt
+from netCDF4 import Dataset
+import numpy as np
+from mpl_toolkits.basemap import Basemap
+
+from pynemo.tests import bdy_coords as bdc
+
+bdc.process_bdy('/Users/thopri/Projects/PyNEMO/inputs/namelist_remote.bdy',False)
+
+rootgrp = Dataset('/Users/thopri/Projects/PyNEMO/outputs/NNA_R12_bdyT_y1979m11.nc', "r", format="NETCDF4")
+bdy_msk = np.squeeze(rootgrp.variables['bdy_msk'][:]) - 1
+bdy_lon = np.squeeze(rootgrp.variables['nav_lon'][:])
+bdy_lat = np.squeeze(rootgrp.variables['nav_lat'][:])
+rootgrp.close()
+
+rootgrp = Dataset('/Users/thopri/Projects/PyNEMO/outputs/coordinates.bdy.nc', "r", format="NETCDF4")
+bdy_it = np.squeeze(rootgrp.variables['nbit'][:])
+bdy_jt = np.squeeze(rootgrp.variables['nbjt'][:])
+bdy_rt = np.squeeze(rootgrp.variables['nbrt'][:])
+rootgrp.close()
+
+bdy_msk = np.ma.masked_where(bdy_msk < 0, bdy_msk)
+for f in range(len(bdy_it)):
+ bdy_msk[bdy_jt[f,],bdy_it[f,]] = bdy_rt[f,]
+
+# Plot to check output
+fig=plt.figure(figsize=(12,12))
+ax=fig.add_subplot(111)
+
+map = Basemap(llcrnrlat=41.,urcrnrlat=65.,\
+ llcrnrlon=-22.,urcrnrlon=25.,\
+ rsphere=(6378137.00,6356752.3142),\
+ resolution='l',projection='lcc',\
+ lat_1=57.,lon_0=-12.5)
+map.drawcoastlines()
+map.drawcountries()
+map.fillcontinents(color='grey')
+map.drawmeridians(np.arange(-25.,25.,2),labels=[0,0,0,1])
+map.drawparallels(np.arange(40.,66.,2),labels=[1,0,0,0])
+
+cmap = plt.cm.get_cmap('jet',10)
+cmaplist = [cmap(i) for i in range(cmap.N)]
+cmaplist[0] = (.7,.7,.7,1.0)
+cmap = cmap.from_list('Custom cmap', cmaplist, cmap.N)
+
+x1,y1 = map(bdy_lon[:,:],bdy_lat[:,:])
+im = map.pcolormesh(x1, y1, bdy_msk, cmap=cmap, vmin=-0.5, vmax=9.5)
+cb = plt.colorbar(orientation='vertical', shrink=0.75, aspect=30, fraction=0.1,pad=0.05)
+cb.set_label('RimWidth Number')
+cb.set_ticks(np.arange(10))
+th=plt.title(('BDY Points'),fontweight='bold')
+
+plt.show()
\ No newline at end of file
diff --git a/test_scripts/bdy_var_plot.py b/test_scripts/bdy_var_plot.py
new file mode 100644
index 0000000000000000000000000000000000000000..cfc47308ab09ac84b0e1adb4238fc78ebaa2ac01
--- /dev/null
+++ b/test_scripts/bdy_var_plot.py
@@ -0,0 +1,259 @@
+#!/usr/bin/env python2
+# -*- coding: utf-8 -*-
+
+import numpy as np
+import scipy.spatial as sp
+from netCDF4 import Dataset
+import matplotlib.pyplot as plt
+from matplotlib.collections import PatchCollection
+from matplotlib.patches import Polygon
+
+
+def nemo_bdy_order(fname):
+ """
+ Determine the ordering and breaks in BDY files to aid plotting.
+ This function takes the i/j coordinates from BDY files and orders them sequentially
+ making it easier to visualise sections along the open boundary. Breaks in the open
+ boundary are also determined (i.e. where the distance between two points > 2**0.5)
+ Args:
+ fname (str) : filename of BDY file
+ Returns:
+ bdy_ind (dict): re-ordered indices
+ bdy_dst (dict): distance (in model coords) between points
+ bdy_brk (dict): location of the break points in the open boundary
+ """
+
+ # open file pointer and extract data
+ rootgrp = Dataset(fname, "r", format="NETCDF4")
+ nbi = np.squeeze(rootgrp.variables['nbidta'][:, :]) - 1 # subtract 1 for python indexing
+ nbj = np.squeeze(rootgrp.variables['nbjdta'][:, :]) - 1
+ nbr = np.squeeze(rootgrp.variables['nbrdta'][:, :]) - 1
+ lon = np.squeeze(rootgrp.variables['nav_lon'][:, :])
+ lat = np.squeeze(rootgrp.variables['nav_lat'][:, :])
+ rootgrp.close()
+
+ rw = np.amax(nbr) + 1
+ bdy_ind = {}
+ bdy_brk = {}
+ bdy_dst = {}
+ nbdy = []
+
+ for r in range(rw):
+ nbdy.append(np.sum(nbr == r))
+
+ # TODO: deal with domain that spans wrap
+
+ # start with outer rim and work in
+
+ for r in range(rw):
+
+ # set initial constants
+
+ ind = nbr == r
+ nbi_tmp = nbi[ind]
+ nbj_tmp = nbj[ind]
+
+ count = 1
+ id_order = np.zeros((nbdy[r], 1), dtype=int) - 1
+ id_order[0,] = 0
+ flag = False
+ mark = 0
+ source_tree = sp.cKDTree(zip(nbi_tmp, nbj_tmp), balanced_tree=False, compact_nodes=False)
+
+ # order bdy entries
+
+ while count < nbdy[r]:
+
+ nn_dist, nn_id = source_tree.query(zip(nbi_tmp[id_order[count - 1]], nbj_tmp[id_order[count - 1]]),
+ k=3, distance_upper_bound=2.9)
+ if np.sum(id_order == nn_id[0, 1]) == 1: # is the nearest point already in the list?
+ if np.sum(id_order == nn_id[0, 2]) == 1: # is the 2nd nearest point already in the list?
+ if flag == False: # then we've found an end point and we can start the search in earnest!
+
+ flag = True
+ id_order[mark] = id_order[count - 1] # reset values
+ id_order[mark + 1:] = -1 # reset values
+ count = mark + 1 # reset counter
+ else:
+ i = 0 # should this be zero?
+ t = count
+ while count == t:
+ if np.sum(id_order == i) == 1:
+ i += 1
+ else:
+ id_order[count] = i
+ flag = False
+ mark = count
+ count += 1
+ elif nn_id[0, 2] == nbdy[r]:
+ i = 0
+ t = count
+ while count == t:
+ if np.sum(id_order == i) == 1:
+ i += 1
+ else:
+ id_order[count] = i
+ flag = False
+ mark = count
+ count += 1
+ else:
+ id_order[count] = nn_id[0, 2]
+ count += 1
+ else:
+ id_order[count] = nn_id[0, 1]
+ count += 1
+
+ bdy_ind[r] = id_order
+ bdy_dst[r] = np.sqrt((np.diff(np.hstack((nbi_tmp[id_order], nbj_tmp[id_order])), axis=0) ** 2).sum(axis=1))
+ bdy_brk[r], = np.where(bdy_dst[r] > 2 ** 0.5)
+ bdy_brk[r] += 1
+ bdy_brk[r] = np.insert(bdy_brk[r], 0, 0) # insert start point
+ bdy_brk[r] = np.insert(bdy_brk[r], len(bdy_brk[r]), len(id_order)) # insert end point
+ bdy_dst[r] = np.insert(np.cumsum(bdy_dst[r]), 0, 0)
+
+ if r == 2: # change to a valid rw number to get a visual output (outer most rw is zero)
+ f, ax = plt.subplots(nrows=1, ncols=1, figsize=(10, 10))
+ plt.scatter(nbi_tmp[bdy_ind[r][:]], nbj_tmp[bdy_ind[r][:]], s=1, marker='.')
+ for t in np.arange(0, len(bdy_ind[r]), 20):
+ plt.text(nbi_tmp[bdy_ind[r][t]], nbj_tmp[bdy_ind[r][t]], t)
+
+ return bdy_ind, bdy_dst, bdy_brk
+
+
+def plot_bdy(fname, bdy_ind, bdy_dst, bdy_brk, varnam, t, rw):
+ """
+ Determine the ordering and breaks in BDY files to aid plotting.
+ This function takes the i/j coordinates from BDY files and orders them sequentially
+ making it easier to visualise sections along the open boundary. Breaks in the open
+ boundary are also determined (i.e. where the distance between two points > 2**0.5)
+ Args:
+ fname (str) : filename of BDY file
+ Returns:
+ bdy_ind (dict): re-ordered indices
+ bdy_dst (dict): distance (in model coords) between points
+ bdy_brk (dict): location of the break points in the open boundary
+ """
+
+ # need to write in a check that either i or j are single values
+
+ rootgrp = Dataset(fname, "r", format="NETCDF4")
+ var = np.squeeze(rootgrp.variables[varnam][t, :])
+ nbr = np.squeeze(rootgrp.variables['nbrdta'][:, :]) - 1
+ var = var[:, nbr == rw]
+
+ # let us use gdept as the central depth irrespective of whether t, u or v
+
+ try:
+ gdep = np.squeeze(rootgrp.variables['deptht'][:, :, :])
+ except KeyError:
+ try:
+ gdep = np.squeeze(rootgrp.variables['gdept'][:, :, :])
+ except KeyError:
+ try:
+ gdep = np.squeeze(rootgrp.variables['depthu'][:, :, :])
+ except KeyError:
+ try:
+ gdep = np.squeeze(rootgrp.variables['depthv'][:, :, :])
+ except KeyError:
+ print 'depth variable not found'
+
+ rootgrp.close()
+
+ jpk = len(gdep[:, 0])
+ nsc = len(bdy_brk[rw][:]) - 1
+
+ dep = {}
+ dta = {}
+
+ # divide data up into sections and re-order
+
+ for n in range(nsc):
+ dta[n] = np.squeeze(var[:, bdy_ind[rw][bdy_brk[rw][n]:bdy_brk[rw][n + 1]]])
+ dep[n] = np.squeeze(gdep[:, bdy_ind[rw][bdy_brk[rw][n]:bdy_brk[rw][n + 1]]])
+
+ # loop over number of sections and plot
+
+ f, ax = plt.subplots(nrows=1, ncols=1, figsize=(11, 4))
+ ax.plot(dta[0][10, :])
+ ax.set_title('BDY points along section: 1, depth level: 10')
+
+ f, ax = plt.subplots(nrows=nsc, ncols=1, figsize=(14, 10 * nsc))
+
+ for n in range(nsc):
+ print('NSC '+str(n))
+ plt.sca(ax[n])
+
+ # from gdep create some pseudo w points
+
+ gdept = dep[n][:, :]
+ coords = np.arange(0, len(gdept[0, :]))
+ gdepw = np.zeros((len(gdept[:, 0]) + 1, len(gdept[0, :])))
+ for z in range(jpk):
+ gdepw[z + 1, :] = gdept[z, :] + (gdept[z, :] - gdepw[z, :])
+
+ gdepvw = np.zeros((len(gdept[:, 0]) + 1, len(gdept[0, :]) + 1))
+
+ # TODO: put in an adjustment for zps levels
+
+ gdepvw[:, 1:-1] = (gdepw[:, :-1] + gdepw[:, 1:]) / 2
+ gdepvw[:, 0] = gdepvw[:, 1]
+ gdepvw[:, -1] = gdepvw[:, -2]
+
+ # create a pseudo bathymetry from the depth data
+
+ bathy = np.zeros_like(coords)
+ mbath = np.sum(dta[n].mask == 0, axis=0)
+
+ for i in range(len(coords)):
+ bathy[i] = gdepw[mbath[i], i]
+
+ bathy_patch = Polygon(np.vstack((np.hstack((coords[0], coords, coords[-1])),
+ np.hstack((np.amax(bathy[:]), bathy, np.amax(bathy[:]))))).T,
+ closed=True,
+ facecolor=(0.8, 0.8, 0), alpha=0, edgecolor=None)
+
+ # Add patch to axes
+ ax[n].add_patch(bathy_patch)
+ ax[n].set_title('BDY points along section: ' + str(n))
+ patches = []
+ colors = []
+
+ for i in range(len(coords)):
+
+ for k in np.arange(jpk - 2, -1, -1):
+
+ if dta[n][k, i] > -10:
+ x = [coords[i] - 0.5, coords[i], coords[i] + 0.5,
+ coords[i] + 0.5, coords[i], coords[i] - 0.5, coords[i] - 0.5]
+ y = [gdepvw[k + 1, i], gdepw[k + 1, i], gdepvw[k + 1, i + 1],
+ gdepvw[k, i + 1], gdepw[k, i], gdepvw[k, i], gdepvw[k + 1, i]]
+
+ polygon = Polygon(np.vstack((x, y)).T, True)
+ patches.append(polygon)
+ colors = np.append(colors, dta[n][k, i])
+
+ # for i in range(len(coords)):
+ # #print(i)
+ # for k in np.arange(jpk - 2, -1, -1):
+ # x = [coords[i] - 0.5, coords[i], coords[i] + 0.5,
+ # coords[i] + 0.5, coords[i], coords[i] - 0.5, coords[i] - 0.5]
+ # y = [gdepvw[k + 1, i], gdepw[k + 1, i], gdepvw[k + 1, i + 1],
+ # gdepvw[k, i + 1], gdepw[k, i], gdepvw[k, i], gdepvw[k + 1, i]]
+ # plt.plot(x, y, 'k-', linewidth=0.1)
+ # plt.plot(coords[i], gdept[k, i], 'k.', markersize=1)
+
+ plt.plot(coords, bathy, '-', color=(0.4, 0, 0))
+ p = PatchCollection(patches, alpha=0.8, edgecolor='none')
+ p.set_array(np.array(colors))
+ ax[n].add_collection(p)
+ f.colorbar(p, ax=ax[n])
+ ax[n].set_ylim((0, np.max(bathy)))
+ ax[n].invert_yaxis()
+
+ return f
+
+fname = '/Users/thopri/Projects/PyNEMO/outputs/NNA_R12_bdyT_y1979m11.nc'
+ind, dst, brk = nemo_bdy_order(fname)
+f = plot_bdy(fname, ind, dst, brk, 'votemper', 0, 0)
+
+plt.show()
\ No newline at end of file