Changeset 5826 for qesdi


Ignore:
Timestamp:
09/10/09 14:19:37 (10 years ago)
Author:
pnorton
Message:

Improved geoplot's behaviour when dealing with variables with axis in the order of lon/lat instead of lat/lon.

Location:
qesdi/geoplot/trunk/lib/geoplot
Files:
1 added
9 edited

Legend:

Unmodified
Added
Removed
  • qesdi/geoplot/trunk/lib/geoplot/array_util.py

    r5403 r5826  
    1212def asColumnVector(a): 
    1313    return N.reshape(a, (a.shape[0], 1)) 
     14 
    1415def asRowVector(a): 
    1516    return N.reshape(a, (1, a.shape[0])) 
  • qesdi/geoplot/trunk/lib/geoplot/grid_builder_base.py

    r5710 r5826  
    8282        gridValues = self._buildGridValues(reducedVar) 
    8383                 
    84         log.debug("After resize:") 
    85         log.debug("y midpoints min =[" + str(gridMidpointY.min()) + \ 
    86                 "] max =[" + str(gridMidpointY.max()) + "]") 
    87         log.debug("x midpoints min =[" + str(gridMidpointX.min()) + \ 
    88                 "] max =[" + str(gridMidpointX.max()) + "]") 
    89          
    90         log.debug("Diff: y [" + str(gridMidpointY.min() - ymid.min()) + \ 
    91                   "][" + str(gridMidpointY.max() - ymid.max()) ) 
    92         log.debug("Diff: x [" + str(gridMidpointX.min() - xmid.min()) + \ 
    93                   "][" + str(gridMidpointX.max() - xmid.max()) +"]") 
    94          
    95         log.debug("self.cdmsVar.shape = %s" % (self.cdmsVar.shape,)) 
    96         log.debug("gridValues.shape = %s" % (gridValues.shape,)) 
    97         log.debug("gridBoundsX.shape = %s" % (gridBoundsX.shape,)) 
     84#        log.debug("After resize:") 
     85#        log.debug("y midpoints min =[" + str(gridMidpointY.min()) + \ 
     86#                "] max =[" + str(gridMidpointY.max()) + "]") 
     87#        log.debug("x midpoints min =[" + str(gridMidpointX.min()) + \ 
     88#                "] max =[" + str(gridMidpointX.max()) + "]") 
     89#         
     90#        log.debug("Diff: y [" + str(gridMidpointY.min() - ymid.min()) + \ 
     91#                  "][" + str(gridMidpointY.max() - ymid.max()) ) 
     92#        log.debug("Diff: x [" + str(gridMidpointX.min() - xmid.min()) + \ 
     93#                  "][" + str(gridMidpointX.max() - xmid.max()) +"]") 
     94#         
     95#        log.debug("self.cdmsVar.shape = %s" % (self.cdmsVar.shape,)) 
     96#        log.debug("gridValues.shape = %s" % (gridValues.shape,)) 
     97#        log.debug("gridBoundsX.shape = %s" % (gridBoundsX.shape,)) 
    9898 
    9999        return Grid(gridBoundsX, gridBoundsY, gridMidpointX, gridMidpointY, gridValues) 
     
    311311        Builds a numpy array of values for each of the grid boxes.  
    312312         
    313         e.g. values[x,y] is the value for grid box x,y. 
     313        e.g. values[y,x] is the value for grid box y, x. 
    314314        """ 
    315315        data = cdmsVar.getValue() 
    316316        missing = cdmsVar.getMissing() 
    317317         
     318        # this data is extracted such that if the axis are in order x/y 
     319        # data[xIndex, yIndex] = value 
     320        # but if they are in y/x it is: 
     321        # data[yIndex, xIndex] = value 
     322        # this function needs to return the value in terms of data[y,x] for the 
     323        # imshow call to work. 
     324 
     325        if not self.__class__._areAxisInOrderYX(cdmsVar): 
     326            data = data.swapaxes(0,1) 
     327                  
    318328        if missing == None: 
    319329            missing = 1e20 
    320          
     330                 
    321331        return MA.masked_values(data, missing) 
    322332         
  • qesdi/geoplot/trunk/lib/geoplot/grid_builder_lat_lon.py

    r5710 r5826  
    3838        xLimits and yLimits must be of the form (low, high) 
    3939        """ 
     40         
     41        #check that the longitude is circular 
     42         
     43        if not self.cdmsVar.getLongitude().isCircular(): 
     44            log.warning("Longitude axis of variable %s isn't circular, this may cause problems subsetting." % (self.cdmsVar.id,)) 
     45            log.warning("self.cdmsVar.getLongitude().isCircular() = %s" % (self.cdmsVar.getLongitude().isCircular(),)) 
     46            log.warning("self.cdmsVar.getLongitude().isCircularAxis() = %s" % (self.cdmsVar.getLongitude().isCircularAxis(),)) 
     47         
    4048#        log.debug('Resizing grid, xLimits: ' + repr(xLimits) + ' yLimits: ' + repr(yLimits)) 
    4149 
     
    4351        latVals = self.cdmsVar.getLatitude()[:] 
    4452         
    45         log.debug("Starting longitude from %s to %s" % ( self.cdmsVar.getLongitude().getValue()[0], self.cdmsVar.getLongitude().getValue()[-1])) 
    46         log.debug("Starting latitude  from %s to %s" % ( self.cdmsVar.getLatitude().getValue()[0] , self.cdmsVar.getLatitude().getValue()[-1])) 
     53#        log.debug("Starting longitude from %s to %s" % ( self.cdmsVar.getLongitude().getValue()[0], self.cdmsVar.getLongitude().getValue()[-1])) 
     54#        log.debug("Starting latitude  from %s to %s" % ( self.cdmsVar.getLatitude().getValue()[0] , self.cdmsVar.getLatitude().getValue()[-1])) 
    4755         
    4856        lonLowHigh = lonVals[0] <= lonVals[-1] 
     
    5260        # is drawn right to the edge. 
    5361        indicator = 'cce' 
    54          
    55         additional_halo = 1 
    56          
     62                 
    5763        if lonLowHigh: 
    5864            lonSelection = (xLimits[0], xLimits[1], indicator) 
     
    125131        return (gridMidpointX, gridMidpointY) 
    126132 
    127     def _isAxisOrderLonLat(self, cdmsVar): 
     133    @staticmethod 
     134    def _areAxisInOrderYX(cdmsVar): 
     135        "Check if the axis on a given variable are longitude/latitude or not" 
     136 
    128137        lon = cdmsVar.getLongitude() 
    129138        lat = cdmsVar.getLatitude() 
    130139         
    131140        if cdmsVar.getAxisIndex(lon) < cdmsVar.getAxisIndex(lat): 
    132             orderIsLonLat = True 
     141            orderIsLatLon = False 
    133142        else: 
    134             orderIsLonLat = False 
     143            orderIsLatLon = True 
    135144         
    136         return orderIsLonLat 
     145        return orderIsLatLon 
    137146 
    138147if __name__ == '__main__': 
  • qesdi/geoplot/trunk/lib/geoplot/grid_builder_national.py

    r5710 r5826  
    171171        return (xAxis, yAxis) 
    172172 
    173  
     173    @staticmethod  
     174    def _areAxisInOrderYX(cdmsVar): 
     175        #!TODO: implement this method 
     176        return True 
    174177    
    175178if __name__ == '__main__': 
  • qesdi/geoplot/trunk/lib/geoplot/grid_builder_rotated.py

    r5710 r5826  
    315315        return (lonrot, latrot) 
    316316 
    317          
    318  
     317    @staticmethod  
     318    def _areAxisInOrderYX(cdmsVar): 
     319        #!TODO: implement this method 
     320        return True 
    319321     
    320322if __name__ == '__main__': 
  • qesdi/geoplot/trunk/lib/geoplot/grid_drawer.py

    r5710 r5826  
    156156            st = time.time() 
    157157             
    158             log.debug("grid_mapUnits.boundsX.shape = %s" % (grid_mapUnits.boundsX.shape,)) 
    159             log.debug("grid_mapUnits.boundsY.shape = %s" % (grid_mapUnits.boundsY.shape,)) 
    160             log.debug("grid_mapUnits.values.shape = %s" % (grid_mapUnits.values.shape,)) 
     158             
     159            log.debug("grid_mapUnits: boundsX.shape = %s, boundsY.shape = %s, values.shape = %s"  
     160                      % (grid_mapUnits.boundsX.shape, grid_mapUnits.boundsY.shape, grid_mapUnits.values.shape)) 
     161             
     162            assert grid_mapUnits.boundsX.shape == grid_mapUnits.boundsY.shape 
    161163             
    162164            values = grid_mapUnits.values 
    163              
    164             assert grid_mapUnits.boundsX.shape == grid_mapUnits.boundsY.shape 
    165              
    166165            boundsShape = grid_mapUnits.boundsX.shape 
    167166             
     
    169168            if values.shape[0] != boundsShape[0] -1 and values.shape[1] != boundsShape[1] -1: 
    170169                 
     170                message = "Value array shape doesn't match the bounds shape."  
     171                 
    171172                #check if the axis were in the opposite order 
    172173                if values.shape[1] == boundsShape[0] -1 and values.shape[0] == boundsShape[1] -1: 
    173                     values = N.transpose(values)            
     174                    message += " It looks like the data may have been transposed." 
     175                    #log.warning("Transposing data to try and fit bounds to values") 
     176                    #values = N.transpose(values) 
     177                 
     178                raise Exception(message) 
     179                                
    174180 
    175181 
     
    213219        yEnd = Y[-1,-1] 
    214220         
    215         log.debug("yStart = %s" % (yStart,)) 
    216         log.debug("yEnd = %s" % (yEnd,)) 
     221        log.debug("yStart = %s, yEnd = %s" % (yStart, yEnd)) 
    217222         
    218223        if yStart < yEnd: 
     
    220225        else: 
    221226            extent=(X.min(), X.max(),Y.min(), Y.max()) 
    222  
     227         
    223228        kwargs.pop('edgecolors') 
    224229        kwargs.pop('linestyle') 
     
    228233         
    229234        log.debug("extent = %s" % (extent,)) 
    230         log.debug("Z.shape = %s" % (Z.shape,)) 
    231235         
    232236        im = axes.imshow(Z, extent=extent, 
  • qesdi/geoplot/trunk/lib/geoplot/layer_drawer.py

    r5760 r5826  
    9595                     frameon=(not self.transparent)) 
    9696         
    97         log.debug("fig.frameon = %s" % (fig.frameon,)) 
    98         log.debug("fig.get_facecolor() = %s" % (fig.get_facecolor(),)) 
    99         log.debug("self.bgcolour = %s" % (self.bgcolour,)) 
    100          
    10197        return fig         
    10298     
  • qesdi/geoplot/trunk/lib/geoplot/layer_drawer_grid.py

    r5760 r5826  
    2222        self._gridDrawer = GridDrawer() 
    2323         
    24     
    25          
    2624        self.cdmsVar = cdmsVar 
    2725        self.showGridLines = showGridLines 
     
    4442            map = self._getMap(xLimits, yLimits) 
    4543            basemap = map.basemap         
     44         
     45#        map = self._getMap(xLimits, yLimits) 
     46#        map.basemap.drawcoastlines(ax=axes, xLimits=xLimits, yLimits=yLimits)      
    4647         
    4748        norm = self._cb.getNormalize(grid=grid) 
  • qesdi/geoplot/trunk/lib/geoplot/tests/unit/test_grid_builder_lat_lon.py

    r5683 r5826  
    7474        valArr = self.builder._buildGridValues(self.tempVar) 
    7575         
     76        if not GridBuilderLatLon._areAxisInOrderYX(self.tempVar): 
     77            missingArr = missingArr.swapaxes(0,1) 
     78         
    7679        nose.tools.assert_equal(missingArr.tolist(9999), valArr.tolist(9999)) 
    7780     
    7881     
    79 class Test_GridBuidlerLatLonWithLonLat(Test_GridBuidlerLatLonBase): 
     82class Test_GridBuidlerLatLonWithLatLon(Test_GridBuidlerLatLonBase): 
    8083    __test__ = True 
    8184     
     
    106109         
    107110     
    108 class Test_GridBuidlerLatLonWithLatLon(Test_GridBuidlerLatLonBase): 
     111class Test_GridBuidlerLatLonWithLonLat(Test_GridBuidlerLatLonBase): 
    109112     
    110113    #Make sure it can cope with the axis reversed 
     
    136139        self.tempVar.setAxisList([self.axisX, self.axisY]) 
    137140        self.tempVar.setMissing(5) 
     141         
    138142 
    139143class Test_GridBuidlerLatLonWithLonHightLow(Test_GridBuidlerLatLonBase): 
    140      
    141     #Make sure it can cope with the axis reversed 
    142144     
    143145    __test__ = True 
Note: See TracChangeset for help on using the changeset viewer.