Fix bug in sea ice analysis with oceanStreams

[xylar]

StreamsFile has no close method, so calls have been removed.

seaice_modelvsobs now resorts to an ocean restart file if no sea-ice restart file is available. This is similar to the behavior of seaice_timeseries.

Fixes PEP8 formatting in sea ice analysis.

[xylar]

@milenaveneziani, these changes are required to get configs/lanl/config.20161117.beta0.A_WCYCL1850S.ne30_oEC_ICG.edison to run on wolf. I think these small bugs were introduced in #115 but I didn't test this run before merging. I found the bug in the process of testing #123.

If you can just verify that configs/lanl/config.20161117.beta0.A_WCYCL1850S.ne30_oEC_ICG.edison runs for you on wolf, I think that would be a good enough test to merge this PR.

[milenaveneziani]

@xylar, @vanroekel: I ran the test mentioned above on wolf, using the newest conda environment (python/anaconda-2.7-climate_new). All the resulting plots look OK, but I get these warnings that remind me what you two have been mentioning recently:

Plotting 2-d maps of MLD climatologies...
Reading files /lustre/scratch2/turquoise/milena/ACME/20161117.beta0.A_WCYCL1850S.ne30_oEC_ICG.edison/run/mpaso.hist.am.timeSeriesStatsMonthly.0001-12-01.nc through /lustre/scratch2/turquoise/milena/ACME/20161117.beta0.A_WCYCL1850S.ne30_oEC_ICG.edison/run/mpaso.hist.am.timeSeriesStatsMonthly.0004-01-01.nc
/turquoise/usr/projects/climate/milena/MPAS-git-repositories/MPAS-Analysis/mpas_analysis/shared/climatology/climatology.py:268: RuntimeWarning: invalid value encountered in divide
  np.nansum(daysarray[monthsValue-1, :, :], axis=0))

Plotting 2-d maps of SSS climatologies...
Reading files /lustre/scratch2/turquoise/milena/ACME/20161117.beta0.A_WCYCL1850S.ne30_oEC_ICG.edison/run/mpaso.hist.am.timeSeriesStatsMonthly.0001-12-01.nc through /lustre/scratch2/turquoise/milena/ACME/20161117.beta0.A_WCYCL1850S.ne30_oEC_ICG.edison/run/mpaso.hist.am.timeSeriesStatsMonthly.0004-01-01.nc
/usr/projects/climate/SHARED_CLIMATE/anaconda_envs/default-2.7/lib/python2.7/site-packages/dask/array/numpy_compat.py:45: RuntimeWarning: invalid value encountered in divide
  x = np.divide(x1, x2, out)

Are they related to the fact that we are using a newer numpy, maybe?

[vanroekel]

In my experimentation with packages, these warnings arise in the newer dask package. Since they don't seem to influence output, I have not looked into them.

[xylar]

Thanks @milenaveneziani!

That error seems to be related to @vanroekel's nearest neighbor interpolation with the latest dask. That will go away with #123 but will be replaced with a boatload of ncremap-related warnings instead. So fun times!

[milenaveneziani]

I wander if we can suppress some of the ncremap/ncks warnings. @czender, do you have any suggestion?
Here are the warnings we get:

  Make ice thickness plots...
Started processing at Wed Mar  1 11:09:09 MST 2017.
Running remap script ncremap from directory /turquoise/usr/projects/climate/SHARED_CLIMATE/anaconda_envs/default-2.7/bin
NCO binaries version "4.6.4" from directory /turquoise/usr/projects/climate/SHARED_CLIMATE/anaconda_envs/default-2.7/bin
Input files in or relative to directory /users/milena/MPAS-git-repositories/MPAS-Analysis
Intermediate/temporary files written to directory /tmp
Output files to directory /usr/projects/climate/milena/20170120.beta0.GMPAS-QU240/clim/regridded
Input data shaped in "mpas"-order, will permute with "ncpdq -a Time,nVertLevels,maxEdges,MaxEdges2,nCategories,ONE,nEdges,nCells"
Input assumed to be MPAS-O/I data: will apply renormalized regridding (--rnr=0.0)
Source and destination grids will both be read from supplied map-file
Input #00: /usr/projects/climate/milena/20170120.beta0.GMPAS-QU240/clim/mpas_iceThickness_FM.nc
PDQ(in)  : /usr/projects/climate/milena/20170120.beta0.GMPAS-QU240/clim/mpas_iceThickness_FM.nc
PDQ(out) : /tmp/ncremap_tmp_pdq.nc.pid214504
Map/Wgt  : /usr/projects/climate/milena/20170120.beta0.GMPAS-QU240/mapping/mpasClimatologyMapping.nc
Regridded: /usr/projects/climate/milena/20170120.beta0.GMPAS-QU240/clim/regridded/mpas_iceThickness_FM.nc
ncks: INFO nco_fl_open() reports current extended filetype = NC_FORMATX_NC3 does not equal previous extended filetype = NC_FORMATX_HDF5. This is expected when NCO is instructed to convert filetypes, i.e., to read from one type and write to another. And when NCO generates grids or templates (which are always netCDF3) when the input file is netCDF4. It is also expected when multi-file operators receive files known to be of different types. However, it could also indicate an unexpected change in input dataset type of which the user should be cognizant.
ncks: WARNING NC_DOUBLE version of "_FillValue" attribute for iceVolumeCell fails isfinite(), value is NaN
ncks: WARNING nco_rgr_map() reports variable "iceVolumeCell" has _FillValue attribute that fails isfinite() (value is NaN) so results are unpredictable. HINT: If regridding fails or values seem weird, retry after converting _FillValue to normal number with, e.g., "ncatted -a _FillValue,iceVolumeCell,m,f,1.0e36 in.nc out.nc"
ncks: WARNING NC_DOUBLE version of "_FillValue" attribute for iceVolumeCell fails isfinite(), value is NaN
Cleaning-up intermediate files...
Quick plots of results from last regridded file:
ncview  /usr/projects/climate/milena/20170120.beta0.GMPAS-QU240/clim/regridded/mpas_iceThickness_FM.nc &
panoply /usr/projects/climate/milena/20170120.beta0.GMPAS-QU240/clim/regridded/mpas_iceThickness_FM.nc &
Elapsed time 0m0s
Started processing at Wed Mar  1 11:09:30 MST 2017.
Running remap script ncremap from directory /turquoise/usr/projects/climate/SHARED_CLIMATE/anaconda_envs/default-2.7/bin
NCO binaries version "4.6.4" from directory /turquoise/usr/projects/climate/SHARED_CLIMATE/anaconda_envs/default-2.7/bin
Input files in or relative to directory /users/milena/MPAS-git-repositories/MPAS-Analysis
Intermediate/temporary files written to directory /tmp
Output files to directory /usr/projects/climate/milena/20170120.beta0.GMPAS-QU240/clim/regridded
Input data shaped in "mpas"-order, will permute with "ncpdq -a Time,nVertLevels,maxEdges,MaxEdges2,nCategories,ONE,nEdges,nCells"
Input assumed to be MPAS-O/I data: will apply renormalized regridding (--rnr=0.0)
Source and destination grids will both be read from supplied map-file
Input #00: /usr/projects/climate/milena/20170120.beta0.GMPAS-QU240/clim/mpas_iceThickness_ON.nc
PDQ(in)  : /usr/projects/climate/milena/20170120.beta0.GMPAS-QU240/clim/mpas_iceThickness_ON.nc
PDQ(out) : /tmp/ncremap_tmp_pdq.nc.pid215848
Map/Wgt  : /usr/projects/climate/milena/20170120.beta0.GMPAS-QU240/mapping/mpasClimatologyMapping.nc
Regridded: /usr/projects/climate/milena/20170120.beta0.GMPAS-QU240/clim/regridded/mpas_iceThickness_ON.nc
ncks: INFO nco_fl_open() reports current extended filetype = NC_FORMATX_NC3 does not equal previous extended filetype = NC_FORMATX_HDF5. This is expected when NCO is instructed to convert filetypes, i.e., to read from one type and write to another. And when NCO generates grids or templates (which are always netCDF3) when the input file is netCDF4. It is also expected when multi-file operators receive files known to be of different types. However, it could also indicate an unexpected change in input dataset type of which the user should be cognizant.
ncks: WARNING NC_DOUBLE version of "_FillValue" attribute for iceVolumeCell fails isfinite(), value is NaN
ncks: WARNING nco_rgr_map() reports variable "iceVolumeCell" has _FillValue attribute that fails isfinite() (value is NaN) so results are unpredictable. HINT: If regridding fails or values seem weird, retry after converting _FillValue to normal number with, e.g., "ncatted -a _FillValue,iceVolumeCell,m,f,1.0e36 in.nc out.nc"
ncks: WARNING NC_DOUBLE version of "_FillValue" attribute for iceVolumeCell fails isfinite(), value is NaN
Cleaning-up intermediate files...
Quick plots of results from last regridded file:
ncview  /usr/projects/climate/milena/20170120.beta0.GMPAS-QU240/clim/regridded/mpas_iceThickness_ON.nc &
panoply /usr/projects/climate/milena/20170120.beta0.GMPAS-QU240/clim/regridded/mpas_iceThickness_ON.nc &
Elapsed time 0m0s

[xylar]

@milenaveneziani and @czender, part of the reason for the warnings is that xarray wants to produce files with _FillValue=nan and it seems to be very difficult to counteract this behavior. I'm not sure why they chose to do this. So these warnings come up anytime we use to_netcdf to write an xarray data set and then remap it with ncremap.

[xylar]

By the way, these warnings are completely unrelated to this PR so we should discuss them somewhere else...

[czender]

@xylar and @milenaveneziani
If you would specify exactly which lines of output you wish omitted, I can probably arrange that. Is it only the WARNING lines, or the two redundant WARNINGs for the same issue, and/or are there other lines you want excised?

[milenaveneziani]

Thanks @czender for the response. I would say maybe we would only want the following:

Output files to directory...
Map/Wgt  : ...
Regridded: ...

I tried to look for a silent option, but I don't think I saw it.

[xylar]

@milenaveneziani, a rather easy solution to this would be to have individual error and log files for the subprocess call to ncremap, and only display their contents to the stdout/stderr if there was an error. This is probably not something I'll add as part of my current work but could add later on.

I'm okay with having an option to suppress output from ncremap sometime in the future but I have definitely benefitted from the current, robust output in some of my debugging of the current setup. So @czender, please retain at least the option of such robust output. Sorry for the conflicting requests.

[czender]

I'll add a verbosity level switch by which the user can turn up/down the printed output. And I'll try to ensure that the WARNINGs only get printed once. It's dangerous to turn-off all warnings by default, so this will be an opt-out situation. On the other hand, the current output is a bit wordy and hasn't been pared down from what I found most useful when developing. Will post again when the verbosity option works.

[czender]

The latest version of ncremap (4.6.5-alpha04 and counting...) has a --vrb_lvl option that can considerably quiet the output to what you might want. It still prints WARNINGs that IMHO users should see. Examples of vrb_lvl=0 to 3 below. 0 prints nothing, 3 prints most everything. I am inclined to set the default level to 2 which shows each file transformation stage inputs and outputs.

Please play around with it and give me feedback within the next week if you want changes before I release 4.6.5 (which will presumably become your default after conda updates).

zender@aerosol:$ ncremap --vrb_lvl=0 -P mpas -i ${DATA}/hdf/mpascice.hist.am.timeSeriesStatsMonthly.0251-01-01.nc -m ${DATA}/maps/map_oEC60to30_to_t62_bilin.20160301.nc -o foo.nc
zender@aerosol:$ ncremap --vrb_lvl=1 -P mpas -i ${DATA}/hdf/mpascice.hist.am.timeSeriesStatsMonthly.0251-01-01.nc -m ${DATA}/maps/map_oEC60to30_to_t62_bilin.20160301.nc -o foo.nc
Input #00: /Users/zender/data/hdf/mpascice.hist.am.timeSeriesStatsMonthly.0251-01-01.nc
Map/Wgt : /Users/zender/data/maps/map_oEC60to30_to_t62_bilin.20160301.nc
Regridded: foo.nc
zender@aerosol:$ ncremap --vrb_lvl=2 -P mpas -i ${DATA}/hdf/mpascice.hist.am.timeSeriesStatsMonthly.0251-01-01.nc -m ${DATA}/maps/map_oEC60to30_to_t62_bilin.20160301.nc -o foo.nc
Input #00: /Users/zender/data/hdf/mpascice.hist.am.timeSeriesStatsMonthly.0251-01-01.nc
att(in) : /Users/zender/data/hdf/mpascice.hist.am.timeSeriesStatsMonthly.0251-01-01.nc
att(out) : /var/folders/p6/qrjddnkd0vb7nkv45mzg48pm0000gn/T/ncremap_tmp_att.nc.pid76522
PDQ(in) : /var/folders/p6/qrjddnkd0vb7nkv45mzg48pm0000gn/T/ncremap_tmp_att.nc.pid76522
PDQ(out) : /var/folders/p6/qrjddnkd0vb7nkv45mzg48pm0000gn/T/ncremap_tmp_pdq.nc.pid76522
Map/Wgt : /Users/zender/data/maps/map_oEC60to30_to_t62_bilin.20160301.nc
Regridded: foo.nc
zender@aerosol:$ ncremap --vrb_lvl=3 -P mpas -i ${DATA}/hdf/mpascice.hist.am.timeSeriesStatsMonthly.0251-01-01.nc -m ${DATA}/maps/map_oEC60to30_to_t62_bilin.20160301.nc -o foo.nc
Started processing at Wed Mar 8 09:07:29 PST 2017.
Running remap script ncremap from directory /Users/zender/bin
NCO binaries version "4.6.5-alpha04" from directory /Users/zender/bin
Input files in or relative to directory /Users/zender
Intermediate/temporary files written to directory /var/folders/p6/qrjddnkd0vb7nkv45mzg48pm0000gn/T
Output files to directory .
Input data shaped in "mpas"-order, will permute with "ncpdq -a Time,nVertLevels,maxEdges,MaxEdges2,nCategories,ONE,nEdges,nCells"
Input assumed to be MPAS-O/I data: will apply renormalized regridding (--rnr=0.0)
Input assumed to be barenaked of attributes: will annotate NC_DOUBLE variables with _FillValue prior to regridding
Source and destination grids will both be read from supplied map-file
Input #00: /Users/zender/data/hdf/mpascice.hist.am.timeSeriesStatsMonthly.0251-01-01.nc
att(in) : /Users/zender/data/hdf/mpascice.hist.am.timeSeriesStatsMonthly.0251-01-01.nc
att(out) : /var/folders/p6/qrjddnkd0vb7nkv45mzg48pm0000gn/T/ncremap_tmp_att.nc.pid76590
PDQ(in) : /var/folders/p6/qrjddnkd0vb7nkv45mzg48pm0000gn/T/ncremap_tmp_att.nc.pid76590
PDQ(out) : /var/folders/p6/qrjddnkd0vb7nkv45mzg48pm0000gn/T/ncremap_tmp_pdq.nc.pid76590
Map/Wgt : /Users/zender/data/maps/map_oEC60to30_to_t62_bilin.20160301.nc
Regridded: foo.nc
Cleaning-up intermediate files...

[milenaveneziani]

Thanks @czender. I could play around with this on edison, while testing various MPAS-Analysis PR's: how do I point to 4.6.5-alpha04 there? (we have been using the default module 4.6.2 in the last few months)

[czender]

it's in ~zender/[bin,lib]_edison