Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Issue#231 #232

Merged
merged 4 commits into from
Mar 7, 2018
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
10 changes: 5 additions & 5 deletions Lib/MV2.py
Original file line number Diff line number Diff line change
Expand Up @@ -13,12 +13,12 @@
from numpy.ma import set_fill_value, shape, size, isMA, isMaskedArray, is_mask, isarray # noqa
from numpy.ma import make_mask, mask_or, nomask # noqa
from numpy import sctype2char, get_printoptions, set_printoptions
from .avariable import AbstractVariable, getNumericCompatibility
from .tvariable import TransientVariable, asVariable
from .grid import AbstractRectGrid
from .error import CDMSError
from cdms2.avariable import AbstractVariable, getNumericCompatibility
from cdms2.tvariable import TransientVariable, asVariable
from cdms2.grid import AbstractRectGrid
from cdms2.error import CDMSError
# from numpy.ma import *
from .axis import allclose as axisAllclose, TransientAxis, concatenate as axisConcatenate, take as axisTake
from cdms2.axis import allclose as axisAllclose, TransientAxis, concatenate as axisConcatenate, take as axisTake


create_mask = make_mask_none
Expand Down
12 changes: 6 additions & 6 deletions Lib/axis.py
Original file line number Diff line number Diff line change
Expand Up @@ -1856,20 +1856,20 @@ def __init__(self, data, bounds=None, id=None,
self._data_ = data[:]
else:
self._data_ = numpy.array(data[:])
elif isinstance(data, numpy.ndarray):
if copy == 0:
self._data_ = data
else:
self._data_ = numpy.array(data)
elif isinstance(data, numpy.ma.MaskedArray):
if numpy.ma.getmask(data) is not numpy.ma.nomask:
if numpy.ma.getmask(data).any() is numpy.bool_(True):
raise CDMSError(
'Cannot construct an axis with a missing value.')
data = data.data
if copy == 0:
self._data_ = data
else:
self._data_ = numpy.array(data)
elif isinstance(data, numpy.ndarray):
if copy == 0:
self._data_ = data
else:
self._data_ = numpy.array(data)
elif data is None:
self._data_ = None
else:
Expand Down
26 changes: 14 additions & 12 deletions Lib/mvCdmsRegrid.py
Original file line number Diff line number Diff line change
Expand Up @@ -445,12 +445,7 @@ def __call__(self, srcVar, **args):

# establish the destination data. Initialize to missing values or 0.
dstData = numpy.ones(dstShape, dtype=srcVar.dtype)
if missingValue is not None and \
re.search('conserv', self.regridMethod) is None:
dstData *= missingValue
else:
dstData *= 0.0

dstData *= missingValue
# sometimes the masked values are not set to missing_values,
# sorry for the extra copy
srcData = numpy.array(
Expand All @@ -463,6 +458,19 @@ def __call__(self, srcVar, **args):
missingValue=missingValue,
**args)

# regrid the mask if we have a mask
if numpy.array(srcVar.mask).any() is True:
dstMask = numpy.ones(dstData.shape)
dstMask[:] *= missingValue
self.regridObj.apply(srcVar.mask, dstMask,
rootPe=0,
missingValue=None,
**args)
elif numpy.any(dstData == missingValue):
# if the missing value is present in the destination data, set
# destination mask
dstMask = (dstData > srcVar.max())

# fill in diagnostic data
if 'diag' in args:
self.regridObj.fillInDiagnosticData(diag=args['diag'], rootPe=0)
Expand All @@ -477,11 +485,6 @@ def __call__(self, srcVar, **args):
if isinstance(v, bytes):
attrs[a] = v

# if the missing value is present in the destination data, set
# destination mask
if numpy.any(dstData == missingValue):
dstMask = (dstData == missingValue)

# create the transient variable. Note: it is unclear whether
# we should create the variable on the supplied dstGrid or
# the local grid.
Expand All @@ -492,5 +495,4 @@ def __call__(self, srcVar, **args):
grid=self.dstGrid,
attributes=attrs,
id=srcVar.id + '_CdmsRegrid')

return dstVar
15 changes: 10 additions & 5 deletions Src/Cdunifmodule.c
Original file line number Diff line number Diff line change
Expand Up @@ -2527,8 +2527,8 @@ PyCdunifVariable_ReadAsArray(PyCdunifVariableObject *self,
array = (PyArrayObject *) PyArray_SimpleNew(d, dims, self->type);
}

if (array != NULL && nitems > 0) {
if (self->nd == 0) {
if (nitems > 0) {
if ((self->nd == 0) && (array != NULL)) {
long zero = 0;
int ret;
Py_BEGIN_ALLOW_THREADS
Expand Down Expand Up @@ -2609,11 +2609,16 @@ PyCdunifVariable_ReadAsArray(PyCdunifVariableObject *self,
}
PyMem_Free(value);
} else {
ret = cdvargets(self->file, self->id, start, count, stride,
array->data);
if(array != NULL) {
ret = cdvargets(self->file,
self->id, start,
count, stride,
array->data);
} else {
ret = -1;
}
}
release_Cdunif_lock()

Py_END_ALLOW_THREADS
;
if (ret == -1) {
Expand Down
60 changes: 30 additions & 30 deletions regrid2/Lib/esmf.py
Original file line number Diff line number Diff line change
Expand Up @@ -167,7 +167,7 @@ def __init__(self, shape, coordSys=ESMF.CoordSys.SPH_DEG,
# ESMF grid object
self.grid = None
# number of cells in [z,] y, x on all processors
self.shape = shape
self.shape = shape[::-1]
# number of dimensions
self.ndims = len(self.shape)
# whether or not cell areas were set
Expand All @@ -177,12 +177,13 @@ def __init__(self, shape, coordSys=ESMF.CoordSys.SPH_DEG,
# whether or not cell centered coordinates were set
self.centersSet = False

maxIndex = numpy.array(shape, dtype=numpy.int32)

# assume last 2 dimensions are Y,X
# For esmf reverse to X, Y
maxIndex = numpy.array(shape[::-1], dtype=numpy.int32)

self.centersSet = False
periodic_dim = self.ndims - 1
pole_dim = self.ndims - 2
periodic_dim = 0
pole_dim = 1
if periodicity == 0:
self.grid = ESMF.Grid(max_index=maxIndex, num_peri_dims=0, staggerloc=[staggerloc],
coord_sys=coordSys)
Expand Down Expand Up @@ -220,7 +221,7 @@ def getLocalSlab(self, staggerloc):
"""
lo, hi = self.getLoHiBounds(staggerloc)
return tuple([slice(lo[i], hi[i], None)
for i in range(self.ndims)])
for i in range(self.ndims)])[::-1]

def getLoHiBounds(self, staggerloc):
"""
Expand All @@ -240,7 +241,7 @@ def getCoordShape(self, staggerloc):
@return tuple
"""
lo, hi = self.getLoHiBounds(staggerloc)
return tuple([hi[i] - lo[i] for i in range(self.ndims)])
return tuple([hi[i] - lo[i] for i in range(self.ndims)])[::-1]

def setCoords(self, coords, staggerloc=CENTER, globalIndexing=False):
"""
Expand All @@ -258,15 +259,14 @@ def setCoords(self, coords, staggerloc=CENTER, globalIndexing=False):
hence the dimensions are reversed here.
"""
# allocate space for coordinates, can only add coordinates once

for i in range(self.ndims):
ptr = self.grid.get_coords(coord_dim=i, staggerloc=staggerloc)
if globalIndexing:
slab = self.getLocalSlab(staggerloc)
slab = self.getLocalSlab(staggerloc)[::-1]
# Populate self.grid with coordinates or the bounds as needed
ptr[:] = coords[self.ndims - i - 1][slab]
ptr[:] = numpy.array(coords[self.ndims - i - 1]).T[slab]
else:
ptr[:] = coords[self.ndims - i - 1]
ptr[:] = numpy.array(coords[self.ndims - i - 1]).T[:]

def getCoords(self, dim, staggerloc):
"""
Expand All @@ -275,8 +275,8 @@ def getCoords(self, dim, staggerloc):
@param staggerloc Stagger location
"""
gridPtr = self.grid.get_coords(coord_dim=dim, staggerloc=staggerloc)
shp = self.getCoordShape(staggerloc)
return numpy.reshape(gridPtr, shp)
shp = self.getCoordShape(staggerloc)[::-1]
return numpy.reshape(gridPtr, shp).T

def setCellAreas(self, areas):
"""
Expand All @@ -287,7 +287,7 @@ def setCellAreas(self, areas):
areaPtr = self.grid.get_item(
item=ESMF.GridItem.AREA,
staggerloc=self.staggerloc)
areaPtr[:] = areas.flat
areaPtr[:] = areas.T.flat
self.cellAreasSet = True

def getCellAreas(self):
Expand All @@ -298,7 +298,7 @@ def getCellAreas(self):
areaPtr = self.grid.get_item(
item=ESMF.GridItem.AREA,
staggerloc=self.staggerloc)
return numpy.reshape(areaPtr, self.shape)
return numpy.reshape(areaPtr, self.shape).T
else:
return None

Expand All @@ -312,7 +312,7 @@ def getMask(self, staggerloc=CENTER):
item=ESMF.GridItem.MASK, staggerloc=staggerloc)
except BaseException:
maskPtr = None
return maskPtr
return maskPtr.T

def setMask(self, mask, staggerloc=CENTER):
"""
Expand All @@ -324,8 +324,8 @@ def setMask(self, mask, staggerloc=CENTER):
maskPtr = self.grid.get_item(
item=ESMF.GridItem.MASK,
staggerloc=staggerloc)
slab = self.getLocalSlab(CENTER)
maskPtr[:] = mask[slab]
slab = self.getLocalSlab(CENTER)[::-1]
maskPtr[:] = mask.T[slab]

def __del__(self):
self.grid.destroy()
Expand Down Expand Up @@ -406,9 +406,9 @@ def getData(self, rootPe):
"""
ptr = self.getPointer()
if rootPe is None:
shp = self.grid.getCoordShape(staggerloc=self.staggerloc)
shp = self.grid.getCoordShape(staggerloc=self.staggerloc)[::-1]
# local data, copy
return numpy.reshape(ptr, shp)
return numpy.reshape(ptr, shp).T
else:
# gather the data on rootPe
lo, hi = self.grid.getLoHiBounds(self.staggerloc)
Expand Down Expand Up @@ -438,7 +438,7 @@ def getData(self, rootPe):
i in range(self.grid.ndims)])
# copy
bigData[slab].flat = ptrs[p]
return bigData # Local
return bigData.T # Local

# rootPe is not None and self.pe != rootPe
return None
Expand All @@ -455,10 +455,10 @@ def setLocalData(self, data, staggerloc, globalIndexing=False):
"""
ptr = self.field.data
if globalIndexing:
slab = self.grid.getLocalSlab(staggerloc)
ptr[:] = data[slab]
slab = self.grid.getLocalSlab(staggerloc)[::-1]
ptr[:] = data.T[slab]
else:
ptr[:] = data
ptr[:] = data.T


##########################################################################
Expand Down Expand Up @@ -559,7 +559,7 @@ def getSrcAreas(self, rootPe):
@return numpy array or None if interpolation is not conservative
"""
if self.srcAreaField is not None:
return self.srcAreaField.data
return self.srcAreaField.data.T
return None

def getDstAreas(self, rootPe):
Expand All @@ -570,7 +570,7 @@ def getDstAreas(self, rootPe):
@return numpy array or None if interpolation is not conservative
"""
if self.srcAreaField is not None:
return self.dstAreaField.data
return self.dstAreaField.data.T
return None

def getSrcAreaFractions(self, rootPe):
Expand All @@ -582,7 +582,7 @@ def getSrcAreaFractions(self, rootPe):
"""
if self.srcFracField is not None:
# self.srcFracField.get_area()
return self.srcFracField.data
return self.srcFracField.data.T
return None

def getDstAreaFractions(self, rootPe):
Expand All @@ -594,7 +594,7 @@ def getDstAreaFractions(self, rootPe):
"""
if self.dstFracField is not None:
# self.dstFracField.get_area()
return self.dstFracField.data
return self.dstFracField.data.T
return None

def __call__(self, srcField=None, dstField=None, zero_region=None):
Expand All @@ -613,8 +613,8 @@ def __call__(self, srcField=None, dstField=None, zero_region=None):

# default is keep the masked values intact
zeroregion = ESMF.Region.SELECT
# if self.regridMethod == CONSERVE:
# zeroregion = None # will initalize to zero
if self.regridMethod == CONSERVE:
zeroregion = None # will initalize to zero

self.regridHandle(
srcfield=srcField.field,
Expand Down
14 changes: 9 additions & 5 deletions regrid2/Lib/mvESMFRegrid.py
Original file line number Diff line number Diff line change
Expand Up @@ -313,19 +313,23 @@ def apply(self, srcData, dstData, rootPe, globalIndexing=False, **args):
# self.maskPtr[:] = srcDataMask[:]
# self.computeWeights(**args)

self.srcFld.field.data[:] = srcData
self.dstFld.field.data[:] = dstData
zero_region = ESMF.Region.SELECT
if 'zero_region' in args.keys():
zero_region=args.get('zero_region')

self.srcFld.field.data[:] = srcData.T
self.dstFld.field.data[:] = dstData.T
# regrid

self.regridObj(self.srcFld.field, self.dstFld.field, zero_region=ESMF.Region.SELECT)
self.regridObj(self.srcFld.field, self.dstFld.field, zero_region=zero_region)

# fill in dstData
if rootPe is None and globalIndexing:
# only fill in the relevant portion of the data
slab = self.dstGrid.getLocalSlab(staggerloc=self.staggerloc)
dstData[slab] = self.dstFld.getData(rootPe=rootPe)
else:
tmp = self.dstFld.field.data
tmp = self.dstFld.field.data.T
if tmp is None:
dstData = None
else:
Expand Down Expand Up @@ -471,7 +475,7 @@ def fillInDiagnosticData(self, diag, rootPe):
'srcAreas', 'dstAreas':
if entry in diag:
diag[entry] = eval(
'self.' + oldMethods[entry] + '(rootPe=rootPe)')
'self.' + oldMethods[entry] + '(rootPe=rootPe)').T
diag['regridTool'] = 'esmf'
diag['regridMethod'] = self.regridMethodStr
diag['periodicity'] = self.periodicity
Expand Down
10 changes: 9 additions & 1 deletion regrid2/Lib/mvGenericRegrid.py
Original file line number Diff line number Diff line change
Expand Up @@ -255,7 +255,7 @@ def apply(self, srcData, dstData,

# set field values to zero where missing, we'll add the mask
# contribution later
indata *= (1 - (srcDataMaskFloat == 1))
# indata *= (1 - (srcDataMaskFloat == 1))

# srcDataMaskFloatData = srcDataMaskFloat * numpy.random.rand(srcHorizShape[0],srcHorizShape[1])*100
# interpolate mask
Expand All @@ -277,6 +277,14 @@ def apply(self, srcData, dstData,
globalIndexing=True,
srcDataMask=srcDataMaskFloat, **args)

# import vcs
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

hum... we should clean this up in the next release.

# pp = vcs.init()
# pp.plot(indata)
# pp.interact()
# pp.clear()
# pp.plot(outdata)
# pp.interact()
# pp.clear()
# apply missing value contribution
if missingValue is not None:
# add mask contribution
Expand Down
2 changes: 2 additions & 0 deletions share/test_data_files.txt
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,8 @@ cc1ce8a03dd04903bf24ae62357b721d v_2000.nc
e1cfa78966764bdd0643ee73654a23fe cdtest14.xml
aab6236b40e5701526da14cf7db304d7 cdtest10.xml
74298ac04f312f0a383c0f4a797d0416 cdtest13.xml
793920296d8ee7f6375bd0f63f98ba02 sftlf.nc
a42e05c76039d90a60a0f3288d826c74 ta.nc
c20b5a0f5f6d031436730df16ad42ab6 tas_mo_clim.nc
acf9323a1f057415cb6a2e06f3a6050d tas_ccsr-95a.xml
4e4648ad9855cb71259926c9d4a1d361 tas_ccsr-95a_1979.01-1979.12.nc
Expand Down
Loading