matplotlib basemap अनुमानों पर चित्रण इलिप्स

Question

मैं बेसमेप प्रोजेक्शन पर इलिप्स को आकर्षित करने की कोशिश कर रहा हूं। पॉलीगॉन जैसे सर्कल को आकर्षित करने के लिए tissot फ़ंक्शन Tissot's indicatrix' को आकर्षित करने के लिए उपयोग किया गया है, जैसा कि निम्न उदाहरण दिखाता है।

from mpl_toolkits.basemap import Basemap 

x0, y0 = 35, -50 
R = 5 

m = Basemap(width=8000000,height=7000000, resolution='l',projection='aea', 
    lat_1=-40.,lat_2=-60,lon_0=35,lat_0=-50) 
m.drawcoastlines() 
m.tissot(x0, y0, R, 100, facecolor='g', alpha=0.5) 

हालांकि, मैं प्रपत्र (x-x0)**2/a**2 + (y-y0)**2/2 = 1 में एक अंडाकार की साजिश रचने में दिलचस्पी है।

import pylab 
from matplotlib.patches import Ellipse 

fig = pylab.figure() 
ax = pylab.subplot(1, 1, 1, aspect='equal') 

x0, y0 = 35, -50 
w, h = 10, 5 
e = Ellipse(xy=(x0, y0), width=w, height=h, linewidth=2.0, color='g') 
ax.add_artist(e) 
e.set_clip_box(ax.bbox) 
e.set_alpha(0.7) 
pylab.xlim([20, 50]) 
pylab.ylim([-65, -35]) 

वहाँ एक tissot के समान प्रभाव के साथ एक आधार मानचित्र प्रक्षेपण पर एक अंडाकार प्लॉट करने के लिए एक रास्ता है: दूसरी ओर, एक नियमित रूप से कार्तीय ग्रिड मैं निम्न नमूना कोड का उपयोग कर सकते पर एक अंडाकार आकर्षित करने के लिए?

Answer 1

बेसमेप के tissot फ़ंक्शन के स्रोत कोड का विश्लेषण करने के घंटों के बाद, ellipses और कुछ डिबगिंग के कुछ गुणों को सीखने के बाद, मैं अपनी समस्या के समाधान के साथ आया। मैं इस प्रकार ellipse नामक एक नया समारोह के साथ आधार मानचित्र वर्ग का विस्तार किया है,

from __future__ import division 
import pylab 
import numpy 

from matplotlib.patches import Polygon 
from mpl_toolkits.basemap import pyproj 
from mpl_toolkits.basemap import Basemap 

class Basemap(Basemap): 
    def ellipse(self, x0, y0, a, b, n, ax=None, **kwargs): 
     """ 
     Draws a polygon centered at ``x0, y0``. The polygon approximates an 
     ellipse on the surface of the Earth with semi-major-axis ``a`` and 
     semi-minor axis ``b`` degrees longitude and latitude, made up of 
     ``n`` vertices. 

     For a description of the properties of ellipsis, please refer to [1]. 

     The polygon is based upon code written do plot Tissot's indicatrix 
     found on the matplotlib mailing list at [2]. 

     Extra keyword ``ax`` can be used to override the default axis instance. 

     Other \**kwargs passed on to matplotlib.patches.Polygon 

     RETURNS 
      poly : a maptplotlib.patches.Polygon object. 

     REFERENCES 
      [1] : http://en.wikipedia.org/wiki/Ellipse 


     """ 
     ax = kwargs.pop('ax', None) or self._check_ax() 
     g = pyproj.Geod(a=self.rmajor, b=self.rminor) 
     # Gets forward and back azimuths, plus distances between initial 
     # points (x0, y0) 
     azf, azb, dist = g.inv([x0, x0], [y0, y0], [x0+a, x0], [y0, y0+b]) 
     tsid = dist[0] * dist[1] # a * b 

     # Initializes list of segments, calculates \del azimuth, and goes on 
     # for every vertex 
     seg = [self(x0+a, y0)] 
     AZ = numpy.linspace(azf[0], 360. + azf[0], n) 
     for i, az in enumerate(AZ): 
      # Skips segments along equator (Geod can't handle equatorial arcs). 
      if numpy.allclose(0., y0) and (numpy.allclose(90., az) or 
       numpy.allclose(270., az)): 
       continue 

      # In polar coordinates, with the origin at the center of the 
      # ellipse and with the angular coordinate ``az`` measured from the 
      # major axis, the ellipse's equation is [1]: 
      # 
      #       a * b 
      # r(az) = ------------------------------------------ 
      #   ((b * cos(az))**2 + (a * sin(az))**2)**0.5 
      # 
      # Azymuth angle in radial coordinates and corrected for reference 
      # angle. 
      azr = 2. * numpy.pi/360. * (az + 90.) 
      A = dist[0] * numpy.sin(azr) 
      B = dist[1] * numpy.cos(azr) 
      r = tsid/(B**2. + A**2.)**0.5 
      lon, lat, azb = g.fwd(x0, y0, az, r) 
      x, y = self(lon, lat) 

      # Add segment if it is in the map projection region. 
      if x < 1e20 and y < 1e20: 
       seg.append((x, y)) 

     poly = Polygon(seg, **kwargs) 
     ax.add_patch(poly) 

     # Set axes limits to fit map region. 
     self.set_axes_limits(ax=ax) 

     return poly 

यह नया समारोह इस उदाहरण की तरह तुरंत इस्तेमाल किया जा सकता:

pylab.close('all') 
pylab.ion() 

m = Basemap(width=12000000, height=8000000, resolution='l', projection='stere', 
      lat_ts=50, lat_0=50, lon_0=-107.) 
m.drawcoastlines() 
m.fillcontinents(color='coral',lake_color='aqua') 
# draw parallels and meridians. 
m.drawparallels(numpy.arange(-80.,81.,20.)) 
m.drawmeridians(numpy.arange(-180.,181.,20.)) 
m.drawmapboundary(fill_color='aqua') 
# draw ellipses 
ax = pylab.gca() 
for y in numpy.linspace(m.ymax/20, 19*m.ymax/20, 9): 
    for x in numpy.linspace(m.xmax/20, 19*m.xmax/20, 12): 
     lon, lat = m(x, y, inverse=True) 
     poly = m.ellipse(lon, lat, 3, 1.5, 100, facecolor='green', zorder=10, 
      alpha=0.5) 
pylab.title("Ellipses on stereographic projection") 

निम्नलिखित में से कौन परिणाम है: