Julia - 在兰伯特等角投影中转换纬度、经度
Posted
技术标签:
【中文标题】Julia - 在兰伯特等角投影中转换纬度、经度【英文标题】:Julia - converting Latitude, Longitude in Lambert conformal projection 【发布时间】:2022-01-23 11:47:57 【问题描述】:我已从 HRRR 下载数据,类似于此笔记本中的 grib2 文件:
https://nbviewer.org/github/microsoft/AIforEarthDataSets/blob/main/data/noaa-hrrr.ipynb
我现在希望将数据用于特定的经度、纬度。但是我不知道如何将我的(经度,纬度)转换为矩阵数据中的网格坐标。
笔记本提到“HRRR 数据来自 Lambert 保形投影”。 (见单元格 8)。
我查看了 GMT 包,他们似乎可以处理 Lambert 保形投影:https://docs.juliahub.com/GMT/EoU0j/0.30.1/proj_examples/。
但是如何转换坐标呢?
下面的代码好像可以转换,但是我觉得这不是给Lambert的,看了GMT documentation之后,无法调整命令中的设置。
lat=37.0; lon=-119.0;
gmt("mapproject -J+proj=merc", [lat;lon])
VectorGMTdataset with 1 segments
First segment DATA
Global BoundingBox: [-1.3247019404399555e7, 4.118821159351122e6]
First seg BoundingBox: [-1.3247019404399555e7, 4.118821159351122e6]
2×1 MatrixFloat64:
4.118821159351122e6
-1.3247019404399555e7
【问题讨论】:
【参考方案1】:发现经纬度其实在grib文件里,所以不用转换:
using GRIB
f = GribFile(grb2_filename)
lons, lats, values = data(Message(f))
# lons in range [225.90452026573686, 299.0828072281622] = [-134.09547973426314, -60.91719277183779]
# lats in range [21.138123000000018, 52.61565330680793]
所以我们可以只查找最近的经纬度索引,读取values中对应的值。
由于纬度和经度都接近110公里,我将距离最小化如下:
(min_error, coord) = findmin(abs.(lats .- lat) .+ abs.(lons .- 360 .- lon))
(0.020456228700048484, CartesianIndex(269, 548))
values[coord]
294.8936767578125
虽然这实际上并没有回答标题问题,但它满足了我当前的需求,也许对其他人有用。
【讨论】:
以上是关于Julia - 在兰伯特等角投影中转换纬度、经度的主要内容,如果未能解决你的问题,请参考以下文章