#
source:
qesdi/geoplot/trunk/lib/geoplot/map_base.py
@
5403

Subversion URL: http://proj.badc.rl.ac.uk/svn/ndg/qesdi/geoplot/trunk/lib/geoplot/map_base.py@5403
Revision 5403, 8.9 KB checked in by pnorton, 12 years ago (diff) |
---|

Line | |
---|---|

1 | """ |

2 | map_base.py |

3 | ============ |

4 | |

5 | A base map class that holds the common map functionality. A map object knows how |

6 | to create a basemap instance and draw it onto a given axes. |

7 | |

8 | """ |

9 | |

10 | import logging |

11 | |

12 | from matplotlib.collections import LineCollection |

13 | import numpy as N |

14 | |

15 | import time |

16 | |

17 | log = logging.getLogger(__name__) |

18 | log.setLevel(logging.DEBUG) |

19 | |

20 | BASEMAP_RESOLUTION = 'i' |

21 | |

22 | |

23 | class MapBase(object): |

24 | """ |

25 | Represents a type of map that can be drawn onto a matplotlib.axes. |

26 | """ |

27 | |

28 | def __init__(self, xLimits, yLimits, drawCoast, drawRivers, addShapefile, **kwargs ): |

29 | """ |

30 | Constructs the map object. |

31 | |

32 | @param xLimits:a tuple representing the (longitude) limits of the map |

33 | @type xLimits: a tuple of 2 floats, corresponding to (min, max) |

34 | @param yLimits:a tuple representing the (latitude) limits of the map |

35 | @type yLimits: a tuple of 2 floats , corresponding to (min, max) |

36 | @param drawCoast: indicates if coaslines should be drawn |

37 | @type drawCoast: boolean |

38 | @param drawRivers: indicates if rivers should be drawn |

39 | @type drawRivers: boolean |

40 | @param addShapefile: location of shapefile to overlay |

41 | @type addShapefile: string |

42 | """ |

43 | |

44 | self.addShapefile = addShapefile |

45 | |

46 | self._drawCoast = drawCoast |

47 | self._drawRivers = drawRivers |

48 | |

49 | self._checkLimits(xLimits, "xLimits") |

50 | self._xLimits = xLimits |

51 | self._checkLimits(yLimits, "yLimits") |

52 | self._yLimits = yLimits |

53 | |

54 | self._setBasemap() |

55 | |

56 | def drawMap(self, axes): |

57 | """ |

58 | Draws the map onto an axes object. |

59 | |

60 | @param axes: the axes on which to draw the map and grid. |

61 | @type axes: an instance of the matplotlib.axes class |

62 | """ |

63 | |

64 | startTime = time.time() |

65 | log.info("Drawing map") |

66 | |

67 | log.debug("limits:" + str(self.xLimits) + ", " + str(self.yLimits)) |

68 | |

69 | log.debug("Values of drawCoast, drawRivers, addShapefile are: " + str(self.drawCoast) + "," + str(self.drawRivers) + "," + str(self.addShapefile)) |

70 | |

71 | # self.basemap.drawmeridians(N.arange(0.,420.,10.), |

72 | # ax = axes, |

73 | # labels=[0,0,0,1]) |

74 | # |

75 | # self.basemap.drawparallels(N.arange(-90.,120.,10.), |

76 | # ax = axes, |

77 | # labels=[1,0,0,0]) |

78 | |

79 | if self.drawCoast: |

80 | |

81 | st = time.time() |

82 | self.basemap.drawcoastlines(ax = axes, linewidth=0.5, color="black", |

83 | antialiased=0, |

84 | xLimits=self.xLimits, yLimits=self.yLimits) |

85 | |

86 | log.debug(" drawn coastlines in %.2fs" % (time.time() - st,)) |

87 | |

88 | if self.drawRivers: |

89 | self.basemap.drawrivers(ax = axes, color='b', linewidth=0.3) |

90 | |

91 | if self.addShapefile != False: |

92 | log.debug("Adding shapefile: " + str(self.addShapefile)) |

93 | self._drawShapefile(axes) |

94 | |

95 | self._configAxes(axes) |

96 | |

97 | self._resizePlot(axes) |

98 | |

99 | log.debug(" drawn map in %.2fs" % (time.time() - startTime,)) |

100 | |

101 | def _setBasemap(self): |

102 | raise NotImplementedError |

103 | |

104 | def _getBasemapResolution(self): |

105 | "Returns the resolution that should be used when creating a basemap" |

106 | |

107 | # if the coast + rivers arn't needed don't load the resolution data, |

108 | # this optimisation is usefull optimisation is basemap takes quite a |

109 | # while to load the coast+river data, and this is skipped if resolution = None |

110 | |

111 | if self.drawCoast == False and self.drawRivers == False: |

112 | resolution = None |

113 | else: |

114 | resolution = BASEMAP_RESOLUTION |

115 | |

116 | log.debug("Using basemap resolution %s" % (resolution,)) |

117 | |

118 | return resolution |

119 | |

120 | def _drawShapefile(self, axes): |

121 | "Uses the self.basemap to draw a shapefile on the axis." |

122 | |

123 | name = 'region_outlines' |

124 | color = 'black' |

125 | linewidth = 0.3 |

126 | |

127 | self.basemap.readshapefile(self.addShapefile, name, |

128 | drawbounds=False) |

129 | |

130 | segments = getattr(self.basemap, name) |

131 | |

132 | reducedSegments = self._cropShapefileSegments(segments) |

133 | |

134 | #add the reduced segments as a new line collection |

135 | col = LineCollection(reducedSegments) |

136 | col.set_color(color) |

137 | col.set_linewidth(linewidth) |

138 | axes.add_collection(col) |

139 | self.basemap.set_axes_limits(ax=axes) |

140 | |

141 | def _getShapefileLimits(self): |

142 | """ |

143 | Returns the lat/lon limits for drawing a shapefile, will not return |

144 | limits larger than SHAPEFILE_LON_LIMITS and SHAPEFILE_LAT_LIMITS but |

145 | if self.xLimits and self.yLimits are smaller than these limits a |

146 | reduced limit will be returned. |

147 | """ |

148 | |

149 | return self.xLimits, self.yLimits |

150 | |

151 | |

152 | def _cropShapefileSegments(self, segments): |

153 | """ |

154 | reduced the shapefile segments to those which have at least 2 |

155 | points within the shapefile limits. |

156 | """ |

157 | lonLimits, latLimits = self._getShapefileLimits() |

158 | |

159 | log.debug("shapefile lonLimits=%s , latLimits = %s" % (lonLimits, latLimits,)) |

160 | |

161 | # get the limits in map projection units |

162 | xLimits, yLimits = self.basemap(lonLimits, latLimits) |

163 | |

164 | log.debug("shapefile xLimits=%s , yLimits = %s" % (xLimits, yLimits,)) |

165 | |

166 | # build a new list of reduced segments |

167 | remainingSegments = [] |

168 | for linePoints in segments: |

169 | |

170 | remainingPoints = [] |

171 | |

172 | for (x, y) in linePoints: |

173 | if self._isPointInLimits(x, y, xLimits, yLimits): |

174 | remainingPoints.append( (x,y) ) |

175 | |

176 | if len(remainingPoints) > 1: |

177 | # To avoid the ends of a line being connected across the map add |

178 | # all the line points not just thoes visible. |

179 | remainingSegments.append(linePoints) |

180 | |

181 | log.debug("number of points remaining = %s" % (sum([len(x) for x in remainingSegments]),)) |

182 | return remainingSegments |

183 | |

184 | def _isPointInLimits(self, x, y, xLimits, yLimits): |

185 | """ |

186 | Checks if a point (x,y) are within the limits xLimits and yLimits, all |

187 | values should be given in map projection units. |

188 | """ |

189 | if x >= xLimits[0] and x <= xLimits[1] and y >=yLimits[0] and y <= yLimits[1]: |

190 | return True |

191 | else: |

192 | return False |

193 | |

194 | def _resizePlot(self, axes): |

195 | """ |

196 | during the drawing process the axes limits may be changed (differently depending |

197 | on what is plotted) this method re-applies the limits of the map so that the expected |

198 | area is displayed in the finished map. |

199 | |

200 | @param axes: the axes on which to draw the map and grid. |

201 | @type axes: an instance of the matplotlib.axes class |

202 | """ |

203 | |

204 | (xMin, yMin) = self.basemap(self.xLimits[0], self.yLimits[0]) |

205 | (xMax, yMax) = self.basemap(self.xLimits[1], self.yLimits[1]) |

206 | |

207 | #log.debug("xMin, xMax:" + str(xMin) + "," + str(xMax)) |

208 | #log.debug("yMin, yMax:" + str(yMin) + "," + str(yMax)) |

209 | |

210 | axes.set_xlim(float(xMin), float(xMax)) |

211 | axes.set_ylim(float(yMin), float(yMax)) |

212 | |

213 | def _checkLimits(self, limits, name): |

214 | |

215 | |

216 | if limits[0] > limits[1]: |

217 | raise Exception("Can't set value of %s for %s, limits should be from low to high"\ |

218 | % (limits, name)) |

219 | |

220 | #getters and setters for the xLimit and yLimit properties, these are needed |

221 | #to try and ensure that the basemap which is dependent on the limits is |

222 | #rebuilt whenever a limit changes. |

223 | |

224 | def __set_xLimits(self,value): |

225 | self._checkLimits(value, "xLimits") |

226 | self._xLimits = value |

227 | self._setBasemap() |

228 | def __get_xLimits(self): return self._xLimits |

229 | xLimits = property(__get_xLimits, __set_xLimits) |

230 | |

231 | def __set_yLimits(self,value): |

232 | self._checkLimits(value, "yLimits") |

233 | self._yLimits = value |

234 | self._setBasemap() |

235 | def __get_yLimits(self): return self._yLimits |

236 | yLimits = property(__get_yLimits, __set_yLimits) |

237 | |

238 | def __set_drawRivers(self, value): |

239 | if value != self._drawRivers: |

240 | self._drawRivers = value |

241 | self._setBasemap() |

242 | |

243 | def __get_drawRivers(self): return self._drawRivers |

244 | drawRivers = property(__get_drawRivers, __set_drawRivers) |

245 | |

246 | def __set_drawCoast(self, value): |

247 | if value != self._drawCoast: |

248 | self._drawCoast = value |

249 | self._setBasemap() |

250 | |

251 | def __get_drawCoast(self): return self._drawCoast |

252 | drawCoast = property(__get_drawCoast, __set_drawCoast) |

253 | |

254 | if __name__ == '__main__': |

255 | #test a simple basemap |

256 | import pylab |

257 | import basemap |

258 | |

259 | bm = basemap.Basemap() |

260 | bm.drawcoastlines() |

261 | pylab.savefig("out.png") |

262 |

**Note:**See TracBrowser for help on using the repository browser.