This is an automated email from the git hooks/post-receive script.

a_valentino-guest pushed a commit to branch master
in repository pyresample.

commit a2064955ae7be4c42373e1942b994a0862c20921
Author: Antonio Valentino <antonio.valent...@tiscali.it>
Date:   Sat Jun 21 10:53:05 2014 +0200

    Imported Upstream version 1.1.0
---
 PKG-INFO                                        |   2 +-
 docs/source/_static/images/uncert_conc_nh.png   | Bin 0 -> 193508 bytes
 docs/source/_static/images/uncert_count_nh.png  | Bin 0 -> 217710 bytes
 docs/source/_static/images/uncert_stddev_nh.png | Bin 0 -> 201293 bytes
 docs/source/conf.py                             |  25 +++
 docs/source/plot.rst                            |   4 +-
 docs/source/swath.rst                           |  37 +++-
 pyresample.egg-info/PKG-INFO                    |   2 +-
 pyresample.egg-info/SOURCES.txt                 |   3 +
 pyresample/kd_tree.py                           | 230 ++++++++++++++++--------
 pyresample/plot.py                              |  42 +----
 pyresample/version.py                           |   2 +-
 test/test_kd_tree.py                            | 144 ++++++++++++---
 test/test_plot.py                               |  18 +-
 14 files changed, 364 insertions(+), 145 deletions(-)

diff --git a/PKG-INFO b/PKG-INFO
index a8e3315..9e9ceba 100644
--- a/PKG-INFO
+++ b/PKG-INFO
@@ -1,6 +1,6 @@
 Metadata-Version: 1.1
 Name: pyresample
-Version: 1.0.0
+Version: 1.1.0
 Summary: Resampling of remote sensing data in Python
 Home-page: UNKNOWN
 Author: Esben S. Nielsen
diff --git a/docs/source/_static/images/uncert_conc_nh.png 
b/docs/source/_static/images/uncert_conc_nh.png
new file mode 100644
index 0000000..56437d4
Binary files /dev/null and b/docs/source/_static/images/uncert_conc_nh.png 
differ
diff --git a/docs/source/_static/images/uncert_count_nh.png 
b/docs/source/_static/images/uncert_count_nh.png
new file mode 100644
index 0000000..90f505f
Binary files /dev/null and b/docs/source/_static/images/uncert_count_nh.png 
differ
diff --git a/docs/source/_static/images/uncert_stddev_nh.png 
b/docs/source/_static/images/uncert_stddev_nh.png
new file mode 100644
index 0000000..0b8f0d5
Binary files /dev/null and b/docs/source/_static/images/uncert_stddev_nh.png 
differ
diff --git a/docs/source/conf.py b/docs/source/conf.py
index d258b7f..8076a20 100644
--- a/docs/source/conf.py
+++ b/docs/source/conf.py
@@ -13,6 +13,31 @@
 
 import sys, os
 
+class Mock(object):
+    def __init__(self, *args, **kwargs):
+        pass
+
+    def __call__(self, *args, **kwargs):
+        return Mock()
+
+    @classmethod
+    def __getattr__(cls, name):
+        if name in ('__file__', '__path__'):
+            return '/dev/null'
+        elif name[0] == name[0].upper():
+            mockType = type(name, (), {})
+            mockType.__module__ = __name__
+            return mockType
+        elif name == "inf":
+            return 0
+        else:
+            return Mock()
+
+MOCK_MODULES = ['numpy', 'pykdtree', 'configobj', 'pyproj', 
+                'scipy', 'scipy.spatial']
+for mod_name in MOCK_MODULES:
+    sys.modules[mod_name] = Mock()
+
 # If extensions (or modules to document with autodoc) are in another directory,
 # add these directories to sys.path here. If the directory is relative to the
 # documentation root, use os.path.abspath to make it absolute, like shown here.
diff --git a/docs/source/plot.rst b/docs/source/plot.rst
index e808d2d..836998d 100644
--- a/docs/source/plot.rst
+++ b/docs/source/plot.rst
@@ -1,3 +1,5 @@
+.. _plot:
+
 Plotting with pyresample and Basemap
 ====================================
 Pyresample supports basic integration with Basemap 
(http://matplotlib.sourceforge.net/basemap).
@@ -133,7 +135,7 @@ Currently only the following set of Proj.4 arguments can be 
interpreted in the c
 {'proj', 'a', 'b', 'ellps', 'lon_0', 'lat_0', 'lon_1', 'lat_1', 'lon_2', 
'lat_2', 'lat_ts'}
 
 Any other Proj.4 parameters will be ignored. 
-If the ellipsoid is not defined in terms of 'ellps', 'a' or ('a', 'b') an 
exception will be raised.
+If the ellipsoid is not defined in terms of 'ellps', 'a' or ('a', 'b') it will 
default to WGS84.
 
 The xsize and ysize in an AreaDefinition will only be used during resampling 
when the image data for use in
 **basemap.imshow** is created. The actual size and shape of the final plot is 
handled by matplotlib.
diff --git a/docs/source/swath.rst b/docs/source/swath.rst
index cfa4bca..8359bce 100644
--- a/docs/source/swath.rst
+++ b/docs/source/swath.rst
@@ -29,10 +29,10 @@ The ImageContainerNearest class can be used for nearest 
neighbour resampling of
  >>> area_con = swath_con.resample(area_def)
  >>> result = area_con.image_data
 
-For other resampling types or splitting the process in two steps use the 
functions in **pyresample.swath** described below. 
+For other resampling types or splitting the process in two steps use the 
functions in **pyresample.kd_tree** described below. 
 
-pyresample.swath
-----------------
+pyresample.kd_tree
+------------------
 
 This module contains several functions for resampling swath data.
 
@@ -164,6 +164,37 @@ If more channels are present in **data** the keyword 
argument **weight_funcs** m
 
 If **data** is a masked array any pixel in the result data that has been 
"contaminated" by weighting of a masked pixel is masked.
 
+Uncertainty estimates
+*********************
+
+Uncertainty estimates in the form of weighted standard deviation can be 
obtained from the **resample_custom** and **resample_gauss** functions.
+By default the functions return the result of the resampling as a single numpy 
array. If the functions are given the keyword argument **with_uncert=True**
+then the following list of numpy arrays will be returned instead: **(result, 
stddev, count)**. **result** is the usual result. **stddev** is the weighted 
standard deviation for each element in the result. **count** is the number of 
data values used in the weighting for each element in the result.
+
+The principle is to view the calculated value for each element in the result 
as a weighted average of values sampled from a statistical variable. 
+An estimate of the standard deviation of the distribution is calculated using 
the unbiased weighted estimator given as 
+**stddev = sqrt((V1 / (V1 ** 2 + V2)) * sum(wi * (xi - result) ** 2))** where 
**result** is the result of the resampling. **xi** is the value of a 
contributing neighbour 
+and **wi** is the corresponding weight. The coefficients are given as **V1 = 
sum(wi)** and **V2 = sum(wi ** 2)**. The standard deviation is only calculated 
for elements in
+the result where more than one neighbour has contributed to the weighting. The 
**count** numpy array can be used for filtering at a higher number of 
contributing neigbours.
+
+Usage only differs in the number of return values from **resample_gauss** and 
**resample_custom**. E.g.:
+
+ >>> result, stddev, count = pr.kd_tree.resample_gauss(swath_def, ice_conc, 
area_def, 
+ ...                                                   
radius_of_influence=20000, 
+ ...                                                   
sigmas=pr.utils.fwhm2sigma(35000), 
+ ...                                                   fill_value=None, 
with_uncert=True)
+
+Below is shown a plot of the result of the resampling using a real data set:
+  .. image:: _static/images/uncert_conc_nh.png
+
+The corresponding standard deviations:
+  .. image:: _static/images/uncert_stddev_nh.png
+
+And the number of contributing neighbours for each element:
+  .. image:: _static/images/uncert_count_nh.png
+
+Notice the standard deviation is only calculated where there are more than one 
contributing neighbour.
+
 Resampling from neighbour info
 ******************************
 The resampling can be split in two steps: 
diff --git a/pyresample.egg-info/PKG-INFO b/pyresample.egg-info/PKG-INFO
index a8e3315..9e9ceba 100644
--- a/pyresample.egg-info/PKG-INFO
+++ b/pyresample.egg-info/PKG-INFO
@@ -1,6 +1,6 @@
 Metadata-Version: 1.1
 Name: pyresample
-Version: 1.0.0
+Version: 1.1.0
 Summary: Resampling of remote sensing data in Python
 Home-page: UNKNOWN
 Author: Esben S. Nielsen
diff --git a/pyresample.egg-info/SOURCES.txt b/pyresample.egg-info/SOURCES.txt
index 8f240ad..99a3d5a 100644
--- a/pyresample.egg-info/SOURCES.txt
+++ b/pyresample.egg-info/SOURCES.txt
@@ -20,6 +20,9 @@ docs/source/_static/images/tb37v_bmng.png
 docs/source/_static/images/tb37v_ortho.png
 docs/source/_static/images/tb37v_pc.png
 docs/source/_static/images/tb37v_quick.png
+docs/source/_static/images/uncert_conc_nh.png
+docs/source/_static/images/uncert_count_nh.png
+docs/source/_static/images/uncert_stddev_nh.png
 pyresample/__init__.py
 pyresample/_multi_proc.py
 pyresample/_spatial_mp.py
diff --git a/pyresample/kd_tree.py b/pyresample/kd_tree.py
index fc1027d..4302a52 100644
--- a/pyresample/kd_tree.py
+++ b/pyresample/kd_tree.py
@@ -89,15 +89,15 @@ def resample_nearest(source_geo_def, data, target_geo_def,
 
 def resample_gauss(source_geo_def, data, target_geo_def,
                    radius_of_influence, sigmas, neighbours=8, epsilon=0,
-                   fill_value=0, reduce_data=True, nprocs=1, segments=None):
+                   fill_value=0, reduce_data=True, nprocs=1, segments=None, 
with_uncert=False):
     """Resamples data using kd-tree gaussian weighting neighbour approach
 
     :Parameters:
     source_geo_def : object
         Geometry definition of source
     data : numpy array               
-        1d array of single channel data points or
-        (source_size, k) array of k channels of datapoints
+        Array of single channel data points or
+        (source_geo_def.shape, k) array of k channels of datapoints
     target_geo_def : object
         Geometry definition of target
     radius_of_influence : float 
@@ -123,14 +123,21 @@ def resample_gauss(source_geo_def, data, target_geo_def,
     segments : {int, None}
         Number of segments to use when resampling.
         If set to None an estimate will be calculated
+    with_uncert : bool, optional
+        Calculate uncertainty estimates
     
-    :Returns: 
-    data : numpy array 
+    :Returns:
+    data : numpy array (default)
         Source data resampled to target geometry
+
+    data, stddev, counts : numpy array, numpy array, numpy array (if 
with_uncert == True)
+        Source data resampled to target geometry.
+        Weighted standard devaition for all pixels having more than one source 
value
+        Counts of number of source values used in weighting per pixel        
     """
-    
+
     def gauss(sigma):
-        #Return gauss functino object
+        #Return gauss function object
         return lambda r: np.exp(-r**2 / float(sigma)**2)
     
     #Build correct sigma argument
@@ -156,20 +163,20 @@ def resample_gauss(source_geo_def, data, target_geo_def,
     return _resample(source_geo_def, data, target_geo_def, 'custom',
                      radius_of_influence, neighbours=neighbours,
                      epsilon=epsilon, weight_funcs=weight_funcs, 
fill_value=fill_value,
-                     reduce_data=reduce_data, nprocs=nprocs, segments=segments)
+                     reduce_data=reduce_data, nprocs=nprocs, 
segments=segments, with_uncert=with_uncert)
 
 def resample_custom(source_geo_def, data, target_geo_def,
                     radius_of_influence, weight_funcs, neighbours=8,
                     epsilon=0, fill_value=0, reduce_data=True, nprocs=1, 
-                    segments=None):
+                    segments=None, with_uncert=False):
     """Resamples data using kd-tree custom radial weighting neighbour approach
 
     :Parameters:
     source_geo_def : object
         Geometry definition of source
     data : numpy array               
-        1d array of single channel data points or
-        (source_size, k) array of k channels of datapoints
+        Array of single channel data points or
+        (source_geo_def.shape, k) array of k channels of datapoints
     target_geo_def : object
         Geometry definition of target
     radius_of_influence : float 
@@ -197,10 +204,16 @@ def resample_custom(source_geo_def, data, target_geo_def,
         Number of segments to use when resampling.
         If set to None an estimate will be calculated
     
-    :Returns: 
-    data : numpy array 
+    :Returns:
+    data : numpy array (default)
         Source data resampled to target geometry
+
+    data, stddev, counts : numpy array, numpy array, numpy array (if 
with_uncert == True)
+        Source data resampled to target geometry.
+        Weighted standard devaition for all pixels having more than one source 
value
+        Counts of number of source values used in weighting per pixel
     """
+
     try:
         for weight_func in weight_funcs:
             if not isinstance(weight_func, types.FunctionType):
@@ -213,11 +226,11 @@ def resample_custom(source_geo_def, data, target_geo_def,
                      radius_of_influence, neighbours=neighbours,
                      epsilon=epsilon, weight_funcs=weight_funcs,
                      fill_value=fill_value, reduce_data=reduce_data,
-                     nprocs=nprocs, segments=segments)
+                     nprocs=nprocs, segments=segments, with_uncert=with_uncert)
 
 def _resample(source_geo_def, data, target_geo_def, resample_type,
              radius_of_influence, neighbours=8, epsilon=0, weight_funcs=None,
-             fill_value=0, reduce_data=True, nprocs=1, segments=None):
+             fill_value=0, reduce_data=True, nprocs=1, segments=None, 
with_uncert=False):
     """Resamples swath using kd-tree approach"""    
                 
     valid_input_index, valid_output_index, index_array, distance_array = \
@@ -230,13 +243,16 @@ def _resample(source_geo_def, data, target_geo_def, 
resample_type,
                                                     nprocs=nprocs,
                                                     segments=segments)
     
+    
     return get_sample_from_neighbour_info(resample_type, 
-                                          target_geo_def.shape, 
-                                          data, valid_input_index, 
-                                          valid_output_index, index_array, 
-                                          distance_array=distance_array, 
-                                          weight_funcs=weight_funcs, 
-                                          fill_value=fill_value)
+                                         target_geo_def.shape, 
+                                         data, valid_input_index, 
+                                         valid_output_index, 
+                                         index_array, 
+                                         distance_array=distance_array, 
+                                         weight_funcs=weight_funcs, 
+                                         fill_value=fill_value, 
+                                         with_uncert=with_uncert)
     
 def get_neighbour_info(source_geo_def, target_geo_def, radius_of_influence, 
                        neighbours=8, epsilon=0, reduce_data=True, nprocs=1, 
segments=None):
@@ -528,11 +544,11 @@ def _create_empty_info(source_geo_def, target_geo_def, 
neighbours):
         
     return valid_output_index, index_array, distance_array 
     
-
 def get_sample_from_neighbour_info(resample_type, output_shape, data, 
                                    valid_input_index, valid_output_index, 
                                    index_array, distance_array=None, 
-                                   weight_funcs=None, fill_value=0):
+                                   weight_funcs=None, fill_value=0, 
+                                   with_uncert=False):
     """Resamples swath based on neighbour info
     
     :Parameters:
@@ -541,6 +557,8 @@ def get_sample_from_neighbour_info(resample_type, 
output_shape, data,
         'custom': Resample based on weight_funcs
     output_shape : (int, int)
         Shape of output as (rows, cols)
+    data : numpy array
+        Source data
     valid_input_index : numpy array
         valid_input_index from get_neighbour_info
     valid_output_index : numpy array
@@ -562,10 +580,10 @@ def get_sample_from_neighbour_info(resample_type, 
output_shape, data,
         with undetermined pixels masked
         
     :Returns: 
-    data : numpy array 
+    result : numpy array 
         Source data resampled to target geometry
     """
-    
+
     if data.ndim > 2 and data.shape[0] * data.shape[1] == 
valid_input_index.size:
         data = data.reshape(data.shape[0] * data.shape[1], data.shape[2])
     elif data.shape[0] != valid_input_index.size:
@@ -575,16 +593,15 @@ def get_sample_from_neighbour_info(resample_type, 
output_shape, data,
         raise ValueError('Mismatch between geometry and dataset')
     
     is_multi_channel = (data.ndim > 1)
-    
     valid_input_size = valid_input_index.sum()
     valid_output_size = valid_output_index.sum()
     
-    if valid_input_size == 0 or valid_output_size == 0:
+    #Handle empty result set
+    if valid_input_size == 0 or valid_output_size == 0: 
         if is_multi_channel:
             output_shape = list(output_shape)
             output_shape.append(data.shape[1])
             
-        #Handle empty result set
         if fill_value is None:
             #Use masked array for fill values
             return np.ma.array(np.zeros(output_shape, data.dtype), 
@@ -592,7 +609,7 @@ def get_sample_from_neighbour_info(resample_type, 
output_shape, data,
         else:
             #Return fill vaues for all pixels
             return np.ones(output_shape, dtype=data.dtype) * fill_value  
-    
+
     #Get size of output and reduced input
     input_size = valid_input_size
     if len(output_shape) > 1:
@@ -644,6 +661,11 @@ def get_sample_from_neighbour_info(resample_type, 
output_shape, data,
         #Add the mask as channels to the dataset
         is_masked_data = True
         new_data = np.column_stack((new_data.data, new_data.mask))
+
+    if new_data.ndim > 1: # Multiple channels or masked input
+        output_shape = list(output_shape)
+        output_shape.append(new_data.shape[1])
+
     
     #Prepare weight_funcs argument for handeling mask data
     if weight_funcs is not None and is_masked_data:
@@ -666,76 +688,131 @@ def get_sample_from_neighbour_info(resample_type, 
output_shape, data,
         result = new_data[new_index_array]
         result[index_mask] = fill_value
     else:
-        #Calculate result using weighting
+        # Calculate result using weighting.
+        # Note: the code below has low readability in order
+        #       to avoid looping over numpy arrays
                 
         #Get neighbours and masks of valid indices
         ch_neighbour_list = []
         index_mask_list = []
-        for i in range(neighbours):
+        for i in range(neighbours): # Iterate over number of neighbours
+            # Make working copy neighbour index and 
+            # set out of bounds indices to zero
             index_ni = index_array[:, i].copy()
             index_mask_ni = (index_ni == input_size)
             index_ni[index_mask_ni] = 0
+
+            # Get channel data for the corresponing indices
             ch_ni = new_data[index_ni]
             ch_neighbour_list.append(ch_ni) 
             index_mask_list.append(index_mask_ni)
         
         #Calculate weights 
         weight_list = []
-        for i in range(neighbours):
-            #Set out of bounds distance to 1 in order to avoid numerical Inf
+        for i in range(neighbours): # Iterate over number of neighbours
+            # Make working copy of neighbour distances and
+            # set out of bounds distance to 1 in order to avoid numerical Inf
             distance = distance_array[:, i].copy()
             distance[index_mask_list[i]] = 1
             
-            if new_data.ndim > 1:
-                #Calculate weights for each channel
-                num_weights = valid_output_index.sum()
+            if new_data.ndim > 1: # More than one channel in data set.
+                # Calculate weights for each channel
                 weights = []
-                for j in range(new_data.shape[1]):                    
+                num_weights = valid_output_index.sum()
+                num_channels = new_data.shape[1]
+                for j in range(num_channels):                    
                     calc_weight = weight_funcs[j](distance)
-                    #Use broadcasting to account for constant weight
+                    # Turn a scalar weight into a numpy array 
+                    # (no effect if calc_weight already is an array)
                     expanded_calc_weight = np.ones(num_weights) * calc_weight
                     weights.append(expanded_calc_weight)
+
+                # Collect weights for all channels for neighbour number    
                 weight_list.append(np.column_stack(weights))
-            else:
+            else: # Only one channel
                 weights = weight_funcs(distance)
                 weight_list.append(weights)
                         
         result = 0
         norm = 0
-        
-        #Calculate result       
-        for i in range(neighbours):   
-            #Find invalid indices to be masked of from calculation
-            if new_data.ndim > 1:
+        count = 0
+        norm_sqr = 0
+        stddev = 0
+
+        # Calculate result       
+        for i in range(neighbours): # Iterate over number of neighbours   
+            # Find invalid indices to be masked of from calculation
+            if new_data.ndim > 1: # More than one channel in data set.
                 inv_index_mask = np.expand_dims(np.invert(index_mask_list[i]), 
axis=1)
-            else:
+            else: # Only one channel
                 inv_index_mask = np.invert(index_mask_list[i])
             
             #Aggregate result and norm
-            result += inv_index_mask * ch_neighbour_list[i] * weight_list[i]
-            norm += inv_index_mask * weight_list[i]
-                                
+            weights_tmp = inv_index_mask * weight_list[i]
+            result += weights_tmp * ch_neighbour_list[i]
+            norm += weights_tmp
+
         #Normalize result and set fillvalue
-        new_valid_index = (norm > 0)
-        result[new_valid_index] /= norm[new_valid_index]
-        result[np.invert(new_valid_index)] = fill_value 
-    
-    #Add fill values
-    if new_data.ndim > 1:
-        full_result = np.ones((output_size, new_data.shape[1])) * fill_value
-    else:
-        full_result = np.ones(output_size) * fill_value
+        result_valid_index = (norm > 0)
+        result[result_valid_index] /= norm[result_valid_index]
+
+        if with_uncert:  # Calculate uncertainties
+            # 2. pass to calculate standard deviation
+            for i in range(neighbours): # Iterate over number of neighbours   
+                # Find invalid indices to be masked of from calculation
+                if new_data.ndim > 1: # More than one channel in data set.
+                    inv_index_mask = 
np.expand_dims(np.invert(index_mask_list[i]), axis=1)
+                else: # Only one channel
+                    inv_index_mask = np.invert(index_mask_list[i])
+
+                # Aggregate stddev information
+                weights_tmp = inv_index_mask * weight_list[i]
+                count += inv_index_mask
+                norm_sqr += weights_tmp ** 2
+                values = inv_index_mask * ch_neighbour_list[i]
+                stddev += weights_tmp * (values - result) ** 2
+
+            # Calculate final stddev
+            new_valid_index = (count > 1)
+            v1 = norm[new_valid_index]
+            v2 = norm_sqr[new_valid_index]
+            stddev[new_valid_index] = np.sqrt((v1 / (v1 ** 2 - v2)) * 
stddev[new_valid_index])
+            stddev[~new_valid_index] = np.NaN
+
+        #Add fill values
+        result[np.invert(result_valid_index)] = fill_value 
+
+    # Create full result
+    if new_data.ndim > 1: # More than one channel
+        output_raw_shape = ((output_size, new_data.shape[1]))
+    else: # One channel
+        output_raw_shape = output_size
+
+    full_result = np.ones(output_raw_shape) * fill_value
     full_result[valid_output_index] = result 
     result = full_result
-    
-    #Calculte correct output shape    
-    if new_data.ndim > 1:
-        output_shape = list(output_shape)
-        output_shape.append(new_data.shape[1])
-    
+
+    if with_uncert:  # Add fill values for uncertainty
+        full_stddev = np.ones(output_raw_shape) * np.nan
+        full_count = np.zeros(output_raw_shape)
+        full_stddev[valid_output_index] = stddev
+        full_count[valid_output_index] = count
+        stddev = full_stddev
+        count = full_count
+        
+        stddev = stddev.reshape(output_shape)
+        count = count.reshape(output_shape)
+
+        if is_masked_data: # Ignore uncert computation of masks
+            stddev = _remask_data(stddev, is_to_be_masked=False)
+            count = _remask_data(count, is_to_be_masked=False)
+
+        # Set masks for invalid stddev
+        stddev = np.ma.array(stddev, mask=np.isnan(stddev))
+
     #Reshape resampled data to correct shape
     result = result.reshape(output_shape)
-    
+
     #Remap mask channels to create masked output
     if is_masked_data:
         result = _remask_data(result)
@@ -746,8 +823,15 @@ def get_sample_from_neighbour_info(resample_type, 
output_shape, data,
         
     #Set output data type to input data type if relevant
     if conserve_input_data_type:
-        result = result.astype(input_data_type)        
-    return result
+        result = result.astype(input_data_type)
+
+    if with_uncert:
+        if np.ma.isMA(result):
+            stddev = np.ma.array(stddev, mask=(result.mask | stddev.mask))
+            count = np.ma.array(count, mask=result.mask)
+        return result, stddev, count
+    else:
+        return result
 
 def _get_fill_mask_value(data_dtype):
     """Returns the maximum value of dtype"""
@@ -761,14 +845,18 @@ def _get_fill_mask_value(data_dtype):
                         data_dtype.type)
     return fill_value
 
-def _remask_data(data):
+def _remask_data(data, is_to_be_masked=True):
     """Interprets half the array as mask for the other half"""
     
     channels = data.shape[-1]
-    mask = data[..., (channels // 2):]            
-    #All pixels affected by masked pixels are masked out
-    mask = (mask != 0)
-    data = np.ma.array(data[..., :(channels // 2)], mask=mask)
+    if is_to_be_masked:
+        mask = data[..., (channels // 2):]            
+        #All pixels affected by masked pixels are masked out
+        mask = (mask != 0)
+        data = np.ma.array(data[..., :(channels // 2)], mask=mask)
+    else:
+        data = data[..., :(channels // 2)]
+
     if data.shape[-1] == 1:
         data = data.reshape(data.shape[:-1])
     return data
diff --git a/pyresample/plot.py b/pyresample/plot.py
index f034015..86e1eed 100644
--- a/pyresample/plot.py
+++ b/pyresample/plot.py
@@ -114,40 +114,14 @@ def area_def2basemap(area_def, **kwargs):
         basemap_args['projection'] = 'cyl'
     else:
         basemap_args['projection'] = area_def.proj_dict['proj']
-    try:
-        basemap_args['lon_0'] = float(area_def.proj_dict['lon_0'])
-    except KeyError:
-        pass
-
-    try:
-        basemap_args['lat_0'] = float(area_def.proj_dict['lat_0']) 
-    except KeyError:
-        pass
-
-    try:
-        basemap_args['lon_1'] = float(area_def.proj_dict['lon_1']) 
-    except KeyError:
-        pass  
-
-    try:
-        basemap_args['lat_1'] = float(area_def.proj_dict['lat_1']) 
-    except KeyError:
-        pass  
-
-    try:
-        basemap_args['lon_2'] = float(area_def.proj_dict['lon_2']) 
-    except KeyError:
-        pass
-
-    try:
-        basemap_args['lat_2'] = float(area_def.proj_dict['lat_2']) 
-    except KeyError:
-        pass
-
-    try:
-        basemap_args['lat_ts'] = float(area_def.proj_dict['lat_ts']) 
-    except KeyError:
-        pass
+    
+    # Try adding potentially remaining args
+    for key in ('lon_0', 'lat_0', 'lon_1', 'lat_1', 'lon_2', 'lat_2', 
+                'lat_ts'):
+        try:
+            basemap_args[key] = float(area_def.proj_dict[key])
+        except KeyError:
+            pass
 
     return Basemap(**basemap_args) 
             
diff --git a/pyresample/version.py b/pyresample/version.py
index c77622e..a278cb1 100644
--- a/pyresample/version.py
+++ b/pyresample/version.py
@@ -15,4 +15,4 @@
 #You should have received a copy of the GNU General Public License
 # along with this program.  If not, see <http://www.gnu.org/licenses/>.
 
-__version__ = '1.0.0'
+__version__ = '1.1.0'
diff --git a/test/test_kd_tree.py b/test/test_kd_tree.py
index 247b5ed..cdd5ef4 100644
--- a/test/test_kd_tree.py
+++ b/test/test_kd_tree.py
@@ -3,14 +3,16 @@ from __future__ import with_statement
 import os
 import sys
 import unittest
+
 import warnings
 if sys.version_info < (2, 6):
     warnings.simplefilter("ignore")
-else:    
+else:
     warnings.simplefilter("always")
+
 import numpy
 
-from pyresample import kd_tree, utils, geometry, grid, data_reduce
+from pyresample import kd_tree, utils, geometry, data_reduce
 
 
 def mp(f):
@@ -45,18 +47,17 @@ class Test(unittest.TestCase):
     tlons = numpy.array([11.280789, 12.649354, 12.080402])
     tlats = numpy.array([56.011037, 55.629675, 55.641535])
     tswath = geometry.SwathDefinition(lons=tlons, lats=tlats)
-    #grid = numpy.ones((1, 1, 2))
-    #grid[0, 0, 0] = 12.562036
-    #grid[0, 0, 1] = 55.715613
     tgrid = geometry.CoordinateDefinition(lons=numpy.array([12.562036]), 
                                           lats=numpy.array([55.715613])) 
-               
+
+                   
     def test_nearest_base(self):     
         res = kd_tree.resample_nearest(self.tswath,\
                                      self.tdata.ravel(), self.tgrid,\
                                      100000, reduce_data=False, segments=1)
         self.assertTrue(res[0] == 2, 'Failed to calculate nearest neighbour')
     
+    @tmp
     def test_gauss_base(self):
         if sys.version_info < (2, 6):
             res = kd_tree.resample_gauss(self.tswath, \
@@ -71,7 +72,8 @@ class Test(unittest.TestCase):
                 self.failIf(('Searching' not in str(w[0].message)), 'Failed to 
create correct neighbour warning')    
         self.assertAlmostEqual(res[0], 2.2020729, 5, \
                                    'Failed to calculate gaussian weighting')
-        
+    
+       
     def test_custom_base(self):
         def wf(dist):
             return 1 - dist/100000.0
@@ -81,7 +83,7 @@ class Test(unittest.TestCase):
                                          self.tdata.ravel(), self.tgrid,\
                                          50000, wf, reduce_data=False, 
segments=1)
         else:
-            with warnings.catch_warnings(record=True) as w:     
+            with warnings.catch_warnings(record=True) as w:
                 res = kd_tree.resample_custom(self.tswath,\
                                              self.tdata.ravel(), self.tgrid,\
                                              50000, wf, reduce_data=False, 
segments=1)
@@ -89,6 +91,56 @@ class Test(unittest.TestCase):
                 self.failIf(('Searching' not in str(w[0].message)), 'Failed to 
create correct neighbour warning')        
         self.assertAlmostEqual(res[0], 2.4356757, 5,\
                                    'Failed to calculate custom weighting')
+
+    @tmp
+    def test_gauss_uncert(self):
+        sigma = utils.fwhm2sigma(41627.730557884883)
+        if sys.version_info < (2, 6):
+            res, stddev, count = kd_tree.resample_gauss(self.tswath, 
self.tdata, 
+                                                         self.tgrid, 100000, 
sigma, 
+                                                         with_uncert=True)
+        else:
+            with warnings.catch_warnings(record=True) as w:
+                res, stddev, count = kd_tree.resample_gauss(self.tswath, 
self.tdata, 
+                                                            self.tgrid, 
100000, sigma, 
+                                                            with_uncert=True)
+                self.failIf(len(w) != 1, 'Failed to create neighbour warning')
+                self.failIf(('Searching' not in str(w[0].message)), 'Failed to 
create correct neighbour warning')
+
+        expected_res = 2.20206560694
+        expected_stddev = 0.707115076173
+        expected_count = 3
+        self.assertAlmostEqual(res[0], expected_res, 5, \
+                                   'Failed to calculate gaussian weighting 
with uncertainty')
+        self.assertAlmostEqual(stddev[0], expected_stddev, 5, \
+                                   'Failed to calculate uncertainty for 
gaussian weighting')
+        self.assertEqual(count[0], expected_count, 'Wrong data point count for 
gaussian weighting with uncertainty')
+
+    @tmp
+    def test_custom_uncert(self):
+        def wf(dist):
+            return 1 - dist/100000.0
+
+        if sys.version_info < (2, 6):
+            res, stddev, counts = kd_tree.resample_custom(self.tswath,
+                                                         self.tdata, 
self.tgrid,
+                                                         100000, wf, 
with_uncert=True)
+        else:
+            with warnings.catch_warnings(record=True) as w:
+                res, stddev, counts = kd_tree.resample_custom(self.tswath,
+                                                         self.tdata, 
self.tgrid,
+                                                         100000, wf, 
with_uncert=True)   
+                self.failIf(len(w) != 1, 'Failed to create neighbour warning')
+                self.failIf(('Searching' not in str(w[0].message)), 'Failed to 
create correct neighbour warning')
+
+        self.assertAlmostEqual(res[0], 2.32193149, 5, \
+                                   'Failed to calculate custom weighting with 
uncertainty')
+        self.assertAlmostEqual(stddev[0], 0.81817972, 5, \
+                                   'Failed to calculate custom for gaussian 
weighting')
+        self.assertEqual(counts[0], 3, 'Wrong data point count for custom 
weighting with uncertainty')
+
+
+    
     def test_nearest(self):
         data = numpy.fromfunction(lambda y, x: y*x, (50, 10))        
         lons = numpy.fromfunction(lambda y, x: 3 + x, (50, 10))
@@ -100,6 +152,7 @@ class Test(unittest.TestCase):
         expected = 15874591.0
         self.assertEqual(cross_sum, expected,\
                              msg='Swath resampling nearest failed')
+    
        
     def test_nearest_1d(self):
         data = numpy.fromfunction(lambda x, y: x * y, (800, 800))        
@@ -115,6 +168,7 @@ class Test(unittest.TestCase):
         self.assertEqual(cross_sum, expected,
                              msg='Swath resampling nearest 1d failed')
     
+    
     def test_nearest_empty(self):
         data = numpy.fromfunction(lambda y, x: y*x, (50, 10))        
         lons = numpy.fromfunction(lambda y, x: 165 + x, (50, 10))
@@ -127,6 +181,7 @@ class Test(unittest.TestCase):
         self.assertEqual(cross_sum, expected,\
                              msg='Swath resampling nearest empty failed')
     
+    
     def test_nearest_empty_multi(self):
         data = numpy.fromfunction(lambda y, x: y*x, (50, 10))        
         lons = numpy.fromfunction(lambda y, x: 165 + x, (50, 10))
@@ -139,6 +194,7 @@ class Test(unittest.TestCase):
         self.assertEqual(res.shape, (800, 800, 3),\
                              msg='Swath resampling nearest empty multi failed')
     
+    
     def test_nearest_empty_multi_masked(self):
         data = numpy.fromfunction(lambda y, x: y*x, (50, 10))        
         lons = numpy.fromfunction(lambda y, x: 165 + x, (50, 10))
@@ -151,7 +207,8 @@ class Test(unittest.TestCase):
                                      fill_value=None)                
         self.assertEqual(res.shape, (800, 800, 3),
                              msg='Swath resampling nearest empty multi masked 
failed')
-            
+    
+    
     def test_nearest_empty_masked(self):
         data = numpy.fromfunction(lambda y, x: y*x, (50, 10))        
         lons = numpy.fromfunction(lambda y, x: 165 + x, (50, 10))
@@ -165,6 +222,7 @@ class Test(unittest.TestCase):
         self.assertTrue(cross_sum == expected,
                         msg='Swath resampling nearest empty masked failed')
     
+    
     def test_nearest_segments(self):
         data = numpy.fromfunction(lambda y, x: y*x, (50, 10))        
         lons = numpy.fromfunction(lambda y, x: 3 + x, (50, 10))
@@ -177,6 +235,7 @@ class Test(unittest.TestCase):
         self.assertEqual(cross_sum, expected,\
                              msg='Swath resampling nearest segments failed')
     
+    
     def test_nearest_remap(self):
         data = numpy.fromfunction(lambda y, x: y*x, (50, 10))        
         lons = numpy.fromfunction(lambda y, x: 3 + x, (50, 10))
@@ -191,7 +250,7 @@ class Test(unittest.TestCase):
         self.assertEqual(cross_sum, expected,\
                              msg='Grid remapping nearest failed')
     
-    @mp
+    
     def test_nearest_mp(self):
         data = numpy.fromfunction(lambda y, x: y*x, (50, 10))        
         lons = numpy.fromfunction(lambda y, x: 3 + x, (50, 10))
@@ -203,7 +262,8 @@ class Test(unittest.TestCase):
         expected = 15874591.0
         self.assertEqual(cross_sum, expected,\
                              msg='Swath resampling mp nearest failed')
-       
+    
+    
     def test_nearest_multi(self):
         data = numpy.fromfunction(lambda y, x: y*x, (50, 10))        
         lons = numpy.fromfunction(lambda y, x: 3 + x, (50, 10))
@@ -217,7 +277,8 @@ class Test(unittest.TestCase):
         expected = 3 * 15874591.0
         self.assertEqual(cross_sum, expected,\
                              msg='Swath multi channel resampling nearest 
failed')
-     
+    
+    
     def test_nearest_multi_unraveled(self):
         data = numpy.fromfunction(lambda y, x: y*x, (50, 10))        
         lons = numpy.fromfunction(lambda y, x: 3 + x, (50, 10))
@@ -230,7 +291,8 @@ class Test(unittest.TestCase):
         expected = 3 * 15874591.0
         self.assertEqual(cross_sum, expected,\
                              msg='Swath multi channel resampling nearest 
failed')
-        
+    
+       
     def test_gauss_sparse(self):
         data = numpy.fromfunction(lambda y, x: y*x, (50, 10))        
         lons = numpy.fromfunction(lambda y, x: 3 + x, (50, 10))
@@ -242,7 +304,8 @@ class Test(unittest.TestCase):
         expected = 15387753.9852
         self.assertAlmostEqual(cross_sum, expected, places=3,\
                                    msg='Swath gauss sparse nearest failed')
-            
+    
+         
     def test_gauss(self):
         data = numpy.fromfunction(lambda y, x: (y + x)*10**-5, (5000, 100))    
    
         lons = numpy.fromfunction(lambda y, x: 3 + (10.0/100)*x, (5000, 100))
@@ -262,7 +325,7 @@ class Test(unittest.TestCase):
         self.assertAlmostEqual(cross_sum, expected,\
                                    msg='Swath resampling gauss failed')
 
-    @tmp
+    
     def test_gauss_fwhm(self):
         data = numpy.fromfunction(lambda y, x: (y + x)*10**-5, (5000, 100))    
    
         lons = numpy.fromfunction(lambda y, x: 3 + (10.0/100)*x, (5000, 100))
@@ -281,7 +344,8 @@ class Test(unittest.TestCase):
         expected = 4872.81050892
         self.assertAlmostEqual(cross_sum, expected,\
                                    msg='Swath resampling gauss failed')
-        
+    
+    
     def test_gauss_multi(self):
         data = numpy.fromfunction(lambda y, x: (y + x)*10**-6, (5000, 100))    
    
         lons = numpy.fromfunction(lambda y, x: 3 + (10.0/100)*x, (5000, 100))
@@ -302,6 +366,40 @@ class Test(unittest.TestCase):
         expected = 1461.84313918
         self.assertAlmostEqual(cross_sum, expected,\
                                    msg='Swath multi channel resampling gauss 
failed')
+
+    @tmp
+    def test_gauss_multi_uncert(self):
+        data = numpy.fromfunction(lambda y, x: (y + x)*10**-6, (5000, 100))    
    
+        lons = numpy.fromfunction(lambda y, x: 3 + (10.0/100)*x, (5000, 100))
+        lats = numpy.fromfunction(lambda y, x: 75 - (50.0/5000)*y, (5000, 100))
+        swath_def = geometry.SwathDefinition(lons=lons, lats=lats)
+        data_multi = numpy.column_stack((data.ravel(), data.ravel(),\
+                                         data.ravel()))
+        if sys.version_info < (2, 6):
+            res, stddev, counts = kd_tree.resample_gauss(swath_def, 
data_multi,\
+                                                self.area_def, 50000, [25000, 
15000, 10000], 
+                                                segments=1, with_uncert=True)
+        else:
+            with warnings.catch_warnings(record=True) as w:
+                res, stddev, counts = kd_tree.resample_gauss(swath_def, 
data_multi,\
+                                                    self.area_def, 50000, 
[25000, 15000, 10000], 
+                                                    segments=1, 
with_uncert=True)
+                self.failIf(len(w) != 1, 'Failed to create neighbour radius 
warning')
+                self.failIf(('Possible more' not in str(w[0].message)), 
'Failed to create correct neighbour radius warning') 
+        cross_sum = res.sum()
+        cross_sum_stddev = stddev.sum()
+        cross_sum_counts = counts.sum()
+        expected = 1461.84313918
+        expected_stddev = 0.446204424799
+        expected_counts = 4934802.0
+        self.assertTrue(res.shape == stddev.shape and stddev.shape == 
counts.shape and counts.shape == (800, 800, 3))
+        self.assertAlmostEqual(cross_sum, expected,
+                                msg='Swath multi channel resampling gauss 
failed on data')
+        self.assertAlmostEqual(cross_sum_stddev, expected_stddev,
+                                msg='Swath multi channel resampling gauss 
failed on stddev')
+        self.assertAlmostEqual(cross_sum_counts, expected_counts,
+                                msg='Swath multi channel resampling gauss 
failed on counts')
+    
     
     def test_gauss_multi_mp(self):
         data = numpy.fromfunction(lambda y, x: (y + x)*10**-6, (5000, 100))    
    
@@ -325,7 +423,8 @@ class Test(unittest.TestCase):
         expected = 1461.84313918
         self.assertAlmostEqual(cross_sum, expected,\
                                    msg='Swath multi channel resampling gauss 
failed') 
-       
+    
+    
     def test_gauss_multi_mp_segments(self):
         data = numpy.fromfunction(lambda y, x: (y + x)*10**-6, (5000, 100))    
    
         lons = numpy.fromfunction(lambda y, x: 3 + (10.0/100)*x, (5000, 100))
@@ -348,7 +447,8 @@ class Test(unittest.TestCase):
         expected = 1461.84313918
         self.assertAlmostEqual(cross_sum, expected,\
                                    msg='Swath multi channel segments 
resampling gauss failed')
-        
+    
+    
     def test_gauss_multi_mp_segments_empty(self):
         data = numpy.fromfunction(lambda y, x: (y + x)*10**-6, (5000, 100))    
    
         lons = numpy.fromfunction(lambda y, x: 165 + (10.0/100)*x, (5000, 100))
@@ -364,6 +464,7 @@ class Test(unittest.TestCase):
                         msg=('Swath multi channel segments empty ' 
                              'resampling gauss failed')) 
     
+    
     def test_custom(self):
         def wf(dist):
             return 1 - dist/100000.0
@@ -385,7 +486,8 @@ class Test(unittest.TestCase):
         expected = 4872.81050729
         self.assertAlmostEqual(cross_sum, expected,\
                                    msg='Swath custom resampling failed')
-     
+    
+    
     def test_custom_multi(self):
         def wf1(dist):
             return 1 - dist/100000.0
@@ -642,6 +744,7 @@ class Test(unittest.TestCase):
         self.assertEqual(cross_sum, expected,\
                              msg='Swath resampling from neighbour info nearest 
failed')
     
+    
     def test_custom_multi_from_sample(self):
         def wf1(dist):
             return 1 - dist/100000.0
@@ -712,9 +815,6 @@ class Test(unittest.TestCase):
         lons = numpy.fromfunction(lambda y, x: 3 + x, (50, 10)) 
         lats = numpy.fromfunction(lambda y, x: 75 - y, (50, 10))
         swath_def = geometry.SwathDefinition(lons=lons, lats=lats)
-#        res = swath.resample_nearest(lons.ravel(), lats.ravel(), 
-#                                    masked_data, self.area_def, 50000,
-#                                    fill_value=None)
         valid_input_index, valid_output_index, index_array, distance_array = \
                                     kd_tree.get_neighbour_info(swath_def, 
                                                              self.area_def, 
diff --git a/test/test_plot.py b/test/test_plot.py
index a1580a5..d513f6e 100644
--- a/test/test_plot.py
+++ b/test/test_plot.py
@@ -5,6 +5,12 @@ import numpy as np
        
 import pyresample as pr
 
+try:
+    import matplotlib
+    matplotlib.use('Agg')
+except ImportError:
+    pass # Postpone fail to individual tests
+
 def tmp(f):
     f.tmp = True
     return f   
@@ -17,8 +23,7 @@ class Test(unittest.TestCase):
     lons = data[:, 0].astype(np.float64)
     lats = data[:, 1].astype(np.float64)
     tb37v = data[:, 2].astype(np.float64)
-
-
+ 
     def test_ellps2axis(self):
         a, b = pr.plot.ellps2axis('WGS84')
         self.assertAlmostEqual(a, 6378137.0, 
@@ -26,7 +31,6 @@ class Test(unittest.TestCase):
         self.assertAlmostEqual(b, 6356752.3142451793, 
                                    msg='Failed to get semi-minor axis of 
ellipsis')
     
-    @tmp   
     def test_area_def2basemap(self):
         area_def = 
pr.utils.parse_area_file(os.path.join(os.path.dirname(__file__), 
                                          'test_files', 'areas.cfg'), 
'ease_sh')[0]
@@ -35,10 +39,7 @@ class Test(unittest.TestCase):
                         bmap.rmajor == 6371228.0, 
                         'Failed to create Basemap object')
 
-               
     def test_plate_carreeplot(self):
-        import matplotlib
-        matplotlib.use('Agg')
         area_def = 
pr.utils.parse_area_file(os.path.join(os.path.dirname(__file__), 
                                             'test_files', 'areas.cfg'), 
'pc_world')[0]
         swath_def = pr.geometry.SwathDefinition(self.lons, self.lats)
@@ -47,10 +48,7 @@ class Test(unittest.TestCase):
                                              fill_value=None)          
         plt = pr.plot._get_quicklook(area_def, result, num_meridians=0, 
                                      num_parallels=0)
-            
     def test_easeplot(self):
-        import matplotlib
-        matplotlib.use('Agg')
         area_def = 
pr.utils.parse_area_file(os.path.join(os.path.dirname(__file__), 
                                             'test_files', 'areas.cfg'), 
'ease_sh')[0]
         swath_def = pr.geometry.SwathDefinition(self.lons, self.lats)
@@ -60,8 +58,6 @@ class Test(unittest.TestCase):
         plt = pr.plot._get_quicklook(area_def, result)
 
     def test_orthoplot(self):
-        import matplotlib
-        matplotlib.use('Agg')
         area_def = 
pr.utils.parse_area_file(os.path.join(os.path.dirname(__file__), 
                                             'test_files', 'areas.cfg'), 
'ortho')[0]
         swath_def = pr.geometry.SwathDefinition(self.lons, self.lats)

-- 
Alioth's /usr/local/bin/git-commit-notice on 
/srv/git.debian.org/git/pkg-grass/pyresample.git

_______________________________________________
Pkg-grass-devel mailing list
Pkg-grass-devel@lists.alioth.debian.org
http://lists.alioth.debian.org/cgi-bin/mailman/listinfo/pkg-grass-devel

Reply via email to