matplotlib / matplotlib/basemap

Discrepancy in LCC (and maybe others) between Lat/Lon corners and Projected Map Boundary

オープン
#646 コメント 0 件 リアクション 0 件 担当者 0 名 GitHub で見る

まだ誰も着手していません。

主要言語
Python
スター
817
フォーク
395
PR マージ指標
30日以内にマージされた PR はありません

説明

### Description:

The reason for this issue is to point out the discrepancy in the behaviour of Basemap when using the Lambert Conformal Conic (LCC) projection with the corner arguments (`llcrnrlon`, `llcrnrlat`, `urcrnrlon`, `urcrnrlat`).

### Issue:

When using an LCC map, I expect the resulting figure to respect the geographic boundaries defined by the Lower-Left (BL) and Upper-Right (TR) corners, e.g. _lines of constant Latitude and of constant Longitude_. However, it appears that Basemap uses the projected X/Y coordinates of these two points and then draws a rectangle in the projected plane connecting them.

Because LCC parallels are curved and meridians converge, a rectangle drawn in the projected X/Y plane does **not** correspond to the geographic box defined by the input coordinates.

The consequence is that data plotted near the boundaries may not be displayed. Parts of the region that are technically within the Latitude/Longitude limits provided are "cut off" because they fall outside the calculated rectangular X/Y box.

### Expected Behaviour:

When a user provides specific lower-left and upper-right coordinates, the expectation is that the resulting map respects those geographic boundaries. I propose two potential solutions:

- **Non-Rectangular Clipping (Cone Section)**: The map boundary (axis spines) should follow the actual meridians and parallels connecting the corners, rather than being forced into a rectangle. This would result in a "cone section" or wedge-shaped plot that accurately reflects the geometry of the projection and the user's requested region.

- **Expand the Projection Box**: Alternatively, if the plot must remain rectangular, the calculated clipping box in the projected plane should be automatically expanded to encompass the maximum extent of the curved geographic region, ensuring that no data strictly within the provided Lat/Lon limits is cut off.

### Example of the issue:

Image

From the image above:

- Red Box: The map boundary created by Basemap (Projected Rectangle).
- Blue Box: The geographic region requested (following Lat/Lon lines).

### MVE:

```python
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
import numpy as np

# --- 1. Map Initialisation ---
fig = plt.figure(figsize=(12, 10))

m = Basemap(
projection='lcc',
resolution='l', # Low resolution
llcrnrlon=-75.0, # Map Lower-Left Longitude
llcrnrlat=25.0, # Map Lower-Left Latitude
urcrnrlon=30.0, # Map Upper-Right Longitude
urcrnrlat=55.0, # Map Upper-Right Latitude
lat_0=55, # Central Latitude (approx mid-point)
lon_0=-30, # Central Longitude (approx mid-point)
lat_1=30, lat_2=60 # Standard Parallels
)

# --- 2. Visual Styling ---
# Pastel colors for clear contrast against the red/blue boxes
water_color = '#AEE1FA'
land_color = '#B2D8B2'

m.drawmapboundary(fill_color=water_color)
m.fillcontinents(color=land_color, lake_color=water_color)
m.drawcoastlines(linewidth=0.6)
m.drawcountries(linewidth=0.5, color='gray')

# Draw grid to visualize curvature
m.drawparallels(np.arange(0, 90, 10), labels=[1, 0, 0, 0], color='gray', linewidth=0.5)
m.drawmeridians(np.arange(-100, 40, 10), labels=[0, 0, 0, 1], color='gray', linewidth=0.5)

# Shared Corner Coordinates
# BL: 40N, 70W | TR: 60N, 0E
bl_lat, bl_lon = 40.0, -70.0
tr_lat, tr_lon = 60.0, 0.0

# Convert Geographic corners to Projected (Meters) corners
bl_x, bl_y = m(bl_lon, bl_lat)
tr_x, tr_y = m(tr_lon, tr_lat)

# --- A. RED BOX: Projected Rectangle (90-degree corners) ---
# This box is defined in the X/Y projection plane.
# It draws straight lines between the projected corners.
# Result: A perfect rectangle on screen, but does not follow lat/lon lines.
proj_x = [bl_x, tr_x, tr_x, bl_x, bl_x]
proj_y = [bl_y, bl_y, tr_y, tr_y, bl_y]

m.plot(proj_x, proj_y, color='red', linewidth=3, label='Projected Box (X/Y Rectangle)')

# --- B. BLUE BOX: Geographic Box (Follows Lat/Lon) ---
# This box is defined by constant Latitude and Longitude.
# Because LCC is conic, these lines curve in the projection.
# We interpolate points to show the true geographic boundary.
lons = np.linspace(bl_lon, tr_lon, 50)
lats = np.linspace(bl_lat, tr_lat, 50)

# Top & Bottom Edges (Constant Latitude)
m.plot(lons, [bl_lat] * 50, latlon=True, color='blue', linewidth=3)
m.plot(lons, [tr_lat] * 50, latlon=True, color='blue', linewidth=3)

# Left & Right Edges (Constant Longitude)
m.plot([bl_lon] * 50, lats, latlon=True, color='blue', linewidth=3)
m.plot([tr_lon] * 50, lats, latlon=True, color='blue', linewidth=3, label='Geographic Box (Lat/Lon)')

# 1. Plot Shared Corners
m.plot(bl_x, bl_y, 'ko', markersize=8, zorder=10)
plt.text(bl_x, bl_y - 150000, 'BL\n(40N, 70W)', ha='center', va='top', fontweight='bold')

m.plot(tr_x, tr_y, 'ko', markersize=8, zorder=10)
plt.text(tr_x, tr_y + 100000, 'TR\n(60N, 0E)', ha='center', va='bottom', fontweight='bold')

# 2. Projected Centre (Red Square)
# The mathematical midpoint in the X/Y plane
mid_x = (bl_x + tr_x) / 2.0
mid_y = (bl_y + tr_y) / 2.0
m.plot(mid_x, mid_y, 'rs', markersize=8, zorder=10)
plt.text(mid_x, mid_y - 200000, 'Proj\nCenter', color='red', ha='center', va='top', fontsize=9, fontweight='bold')

# 3. Geographic Centre (Blue Triangle)
# The mathematical midpoint in Degrees, then projected
mid_lon = (bl_lon + tr_lon) / 2.0
mid_lat = (bl_lat + tr_lat) / 2.0
geo_mid_x, geo_mid_y = m(mid_lon, mid_lat)
m.plot(geo_mid_x, geo_mid_y, 'b^', markersize=8, zorder=10)
plt.text(geo_mid_x, geo_mid_y + 150000, 'Geo\nCenter', color='blue', ha='center', va='bottom', fontsize=9,
fontweight='bold')

plt.legend(loc='lower right')
plt.show()
```

Environment:

- Python version: `3.12.8`
- Matplotlib version: `3.10.0`
- Basemap version: `1.4.1`

コントリビューションガイド

このリポジトリのコントリビューションガイドは索引されていません

はじめの一歩

  1. issue を最後まで読み、次にプロジェクトのコントリビューションガイドを読みます。
  2. 着手することを issue にコメントします — 二人が同じ作業をするのを防げます。
  3. リポジトリをフォークし、ブランチを切って変更します。
  4. issue 番号を参照したプルリクエストを送ります。

調査の方向性

提供された MVE を Basemap の LCC 投影とその隅の引数で再現し、続いて Basemap の初期化と drawmapboundary の動作によって地図の境界がどのように設定されるかを調べます。投影された矩形と補間された緯度/経度の辺を比較します; 合意された境界または拡張の動作によって、有効な地理的領域データがクリップされなくなった場合にのみ、issue は完了とします。

索引モデルが issue の本文から書いたものです。

評価

技術スタック
matplotlib, numpy, python
領域
data-visualization
issue の種類
バグ
難易度
5/5
見積もり時間
1週間以上
活発さ
停滞
明瞭さ
おおむね明確
初心者へのやさしさ
35/100

新しい issue をメールで受け取る

初心者向けの GitHub issue を短くまとめたダイジェスト。