{
"cells": [
{
"cell_type": "markdown",
"id": "100ed2b5",
"metadata": {},
"source": [
"## Parameter `use_crs`\n",
"\n",
"The `use_crs` parameter is supported by all geometry layers which have the `map` parameter, i.e. by layers which can be used for vizualization of Geopandas GeoDataFrame.\n",
"\n",
"Those layers are: `point, path, polygon, text, label and map`.\n",
"\n",
"By default, Lets-Plot converts all coordinates in every GeoDataFrame to degrees of longitude and latitude and leter projects these map coordinates to screen coordinates useing Mercator projection.\n",
"\n",
"Sometimes it might be necessary to use some other coordinate reference system (CRS) instead of the default WGS84 / Mercator.\n",
"\n",
"If this is the case, you can use the `use_crs` parameter to specify EPSG code of CRS to project GeoDataFrame coordinates to."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "b42d401a",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"\n",
" \n",
" \n",
" "
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"import numpy as np\n",
"from lets_plot import *\n",
"LetsPlot.setup_html()"
]
},
{
"cell_type": "markdown",
"id": "785bfaa5",
"metadata": {},
"source": [
"Let's say we have a raster image (for example, a satellite or aerial photograph, or a DEM image) that covers a square lend area with:\n",
"- left-bottom corner at coordinates (25_000, 4_500_000)\n",
"- right-top corner at coordinates (525_000, 5_000_000)\n",
"\n",
"In practice, the image coordinates and the image coordinate reference system (CRS) are known from \n",
"a metadata accompanying the image.\n",
"\n",
"Here, for the sake of simplicity, we will use just a simple generated image 2 x 2 pix.\n",
"\n",
"We will also presume that our image coordinates use [UTM zone 33N](https://epsg.io/32633) CRS (aka EPSG:32633).\n",
"\n",
"Thus, the image coordinates are given in **meters** instead of more familiar longitude and latatude."
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "31d95b93",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
" \n",
" "
],
"text/plain": [
""
]
},
"execution_count": 2,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# image 2 x 2 pix\n",
"image_arr = np.array([\n",
" [[150, 0, 0], [0, 150, 0]],\n",
" [[0, 0, 150], [150, 150, 0]]\n",
"])\n",
"\n",
"image_CRS = \"EPSG:32633\"\n",
"\n",
"# Image coordinates\n",
"min_x = 25_000\n",
"max_x = 525_000\n",
"min_y = 4_500_000\n",
"max_y = 5_000_000\n",
"\n",
"ggplot() + geom_imshow(image_arr, extent=[min_x, max_x, min_y, max_y])"
]
},
{
"cell_type": "markdown",
"id": "06cc7124",
"metadata": {},
"source": [
"### Combining image and map\n",
"\n",
"Step 1 here would be to aquare a Geopandas GeoDataFrame containing map of the area of interesst.\n",
"\n",
"In this demo we will use GeoDataFrame containing a map of Italy which we will obtain using Lets-Plot geocoding module."
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "7ed564da",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"The geodata is provided by © OpenStreetMap contributors and is made available here under the Open Database License (ODbL).\n"
]
}
],
"source": [
"from lets_plot.geo_data import *"
]
},
{
"cell_type": "code",
"execution_count": 4,
"id": "7315f07a",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
" \n",
" "
],
"text/plain": [
""
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"italy = geocode_countries('italy').inc_res().get_boundaries()\n",
"\n",
"ggplot() + geom_map(map=italy)"
]
},
{
"cell_type": "code",
"execution_count": 5,
"id": "b346a82b",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"\n",
"Name: WGS 84\n",
"Axis Info [ellipsoidal]:\n",
"- Lat[north]: Geodetic latitude (degree)\n",
"- Lon[east]: Geodetic longitude (degree)\n",
"Area of Use:\n",
"- name: World.\n",
"- bounds: (-180.0, -90.0, 180.0, 90.0)\n",
"Datum: World Geodetic System 1984 ensemble\n",
"- Ellipsoid: WGS 84\n",
"- Prime Meridian: Greenwich"
]
},
"execution_count": 5,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"italy.crs"
]
},
{
"cell_type": "markdown",
"id": "de48d4db",
"metadata": {},
"source": [
"#### The \"Italy map\" and our image use different CRS-s\n",
"\n",
"The image uses UTM CRS with coordinates in meters and the \"Italy map\" uses the EPSG:4326 CRS with coordinates in degrees of longitude and latitude.\n",
"\n",
"Due to this discrepancy we can not right away combine these two layes in one plot. "
]
},
{
"cell_type": "markdown",
"id": "5c77ee74",
"metadata": {},
"source": [
"#### The `use_crs` parameter\n",
"\n",
"The solution here is actually pretty straightforward: we will just let the map layer know which exactly CRS the other plot layer uses (i.e. which CRS our image uses). \n"
]
},
{
"cell_type": "code",
"execution_count": 6,
"id": "71ec2b02",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
" \n",
" "
],
"text/plain": [
""
]
},
"execution_count": 6,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"ggplot() + geom_map(map=italy, use_crs=image_CRS) + labs(x=\"\", y=\"\")"
]
},
{
"cell_type": "code",
"execution_count": 7,
"id": "1fbb1a33",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
" \n",
" "
],
"text/plain": [
""
]
},
"execution_count": 7,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Image and map\n",
"( ggplot() \n",
" + geom_imshow(image_arr, extent=[min_x, max_x, min_y, max_y])\n",
" + geom_map(map=italy, use_crs=image_CRS) \n",
" + labs(x=\"\", y=\"\")\n",
")"
]
},
{
"cell_type": "code",
"execution_count": 8,
"id": "7fc54352",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
" \n",
" "
],
"text/plain": [
""
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Select better theme and colors\n",
"( ggplot() \n",
" + geom_imshow(image_arr, extent=[min_x, max_x, min_y, max_y])\n",
" + geom_map(map=italy, use_crs=image_CRS, color=\"#a6b6ba\", size=.8) \n",
" + labs(x=\"\", y=\"\")\n",
" + theme_bw()+ flavor_solarized_dark()\n",
" + ggsize(500, 500)\n",
")"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.7.13"
}
},
"nbformat": 4,
"nbformat_minor": 5
}