0

I am focusing on visualizing some datasets of Switzerland through 'contextily' package in python. (Here is just the sample. In fact, the whole dataset has 16000 records in Switzerland.)

Unnamed: 0 User_location_code Twitter_created_time TRUE
16 226333 CH 1/1/2018 (46.2017559, 6.1466014)
23 4488879 CH 1/1/2018 (47.3744489, 8.5410422)
28 3772151 CH 1/1/2018 (47.1486137, 8.5539378)
35 3226920 CH 1/1/2018 (46.4978958, 9.8392428)
51 4577072 CH 1/1/2018 (46.2131843, 6.0815756)
55 4464950 CH 1/1/2018 (46.4602978, 6.8418655)
65 173633 CH 1/1/2018 (46.2767572, 6.168958)
86 4670064 CH 1/1/2018 (46.2922522, 7.5323195)
89 529343 CH 1/1/2018 (46.2017559, 6.1466014)
90 879080 CH 1/2/2018 (47.2617701, 8.4197059)
100 4195443 CH 1/2/2018 (46.2017559, 6.1466014)
106 1837089 CH 1/2/2018 (46.2017559, 6.1466014)
109 4392712 CH 1/2/2018 (46.4705779, 6.8436858)
110 4191284 CH 1/2/2018 (46.6932583, 6.4611912)
112 3581491 CH 1/2/2018 (47.3744489, 8.5410422)
129 4238370 CH 1/3/2018 (46.5332403, 9.874586)
146 3103734 CH 1/4/2018 (46.4602978, 6.8418655)
147 3103733 CH 1/4/2018 (46.4602978, 6.8418655)
179 4531534 CH 1/6/2018 (47.4634767, 7.8125639)
184 207311 CH 1/6/2018 (47.0505452, 8.3054682)
194 974817 CH 1/6/2018 (46.2017559, 6.1466014)
205 4623486 CH 1/7/2018 (46.4978958, 9.8392428)
214 2782686 CH 1/7/2018 (46.9886931, 6.8906893)
221 4394683 CH 1/8/2018 (46.5218269, 6.6327025)
223 586320 CH 1/8/2018 (47.31874, 7.6698284)
229 787444 CH 1/8/2018 (47.4991723, 8.7291498)
238 1615538 CH 1/9/2018 (46.6356963, 6.5320717)
248 543950 CH 1/9/2018 (47.4250593, 9.3765878)
249 1005935 CH 1/9/2018 (47.3744489, 8.5410422)
252 4217738 CH 1/10/2018 (46.1839156, 6.1224047)
269 411809 CH 1/11/2018 (46.2017559, 6.1466014)
271 1661548 CH 1/11/2018 (47.33490645, 9.163183508014553)
279 814333 CH 1/12/2018 (46.6356963, 6.5320717)
286 652801 CH 1/12/2018 (47.3744489, 8.5410422)
293 591620 CH 1/13/2018 (46.2231863, 6.2242767)
297 4378107 CH 1/13/2018 (46.4978958, 9.8392428)
301 2356500 CH 1/13/2018 (46.5939043, 7.9078016)
308 652799 CH 1/13/2018 (47.3744489, 8.5410422)
315 649442 CH 1/14/2018 (46.6583588, 7.6965086)
319 4380071 CH 1/14/2018 (46.2017559, 6.1466014)
320 190831 CH 1/14/2018 (46.5032766, 6.6867879)
321 4384644 CH 1/14/2018 (46.4602978, 6.8418655)
325 4629311 CH 1/14/2018 (46.647529399999996, 7.080368826694413)
339 411808 CH 1/16/2018 (46.2017559, 6.1466014)
341 2783072 CH 1/16/2018 (46.5218269, 6.6327025)
344 4395955 CH 1/16/2018 (46.6789116, 7.1027113)
348 4394190 CH 1/16/2018 (46.2017559, 6.1466014)
354 440279 CH 1/17/2018 (46.7410614, 7.2193275)
358 4397447 CH 1/17/2018 (46.7899498, 7.2183078)

Because of the geographical coordinate in meter, I decided to use the 'EPSG:3857' as the CRS.

# Load the environment
import numpy
import pandas
import geopandas
import pysal
import seaborn
import contextily
import matplotlib.pyplot as plt

# Read the dataset
CH_points = gpd.read_file('CH_points04.shp')

# Print out the crs of the dataset
print(CH_points.crs)

CH_points_3857 = CH_points.to_crs("epsg:3857")

# Print out the crs of the dataset
print(CH_points_3857.crs)

However, I met a mistake during the basemap time. Specifically, in the following code, I have finished the preliminary visualization about the dataset that I used.

# Generate scatter plot for quick checking
seaborn.jointplot(x = "Longitude", y = "Latitude", data = CH_points, s = 2.0);

The result is here, it works well.

And then, I consider using the function 'contextily.add_basemap()' to make kde-plot looks better, at least in Switzerland.

# Set up figure and axis
f, ax = plt.subplots(1, figsize=(10, 10))
seaborn.kdeplot(
    CH_points_3857["X"],
    CH_points_3857["Y"],
    n_levels=50,
    shade=True,
    alpha=0.55,
    cmap="viridis_r",
)
# Add basemap
contextily.add_basemap(
    ax, crs = CH_points_3857.crs.to_string(), source=contextily.providers.OpenStreetMap.CH
)
# Remove axes
ax.set_axis_off()

The result is here. As you see, it can load the basemap, not bad, right? But the basemap now is in France, not in Switzerland.

How can I align these properly?

Note: For pictures, you could see here. https://github.com/geopandas/contextily/issues/212

For the dataset, you could find it here. https://drive.google.com/file/d/1c85ofyiHKDxjyV5m9LgvWNWL_E9fiifb/view?usp=sharing

  • Can you provide the data 'CH_points', 'CH_points_3857'? – Ori Yarden PhD Feb 23 '23 at 18:22
  • Did you mean to include a result? Also, is the data actually in EPSG:3857 or does it need to be transformed? Please include a full [mre] e.g. using sample points so we can help debug. – Michael Delgado Feb 23 '23 at 22:19
  • @OriYarden Thanks a lot. I have edited the question again. And I think you could download the dataset now. – Xiaoyan Lan Feb 26 '23 at 16:48
  • @MichaelDelgado In fact, without enough reputation of stack overflow, I cannot upload the pictures into here. However, thanks for your suggestions, I have edited the description of my question now. – Xiaoyan Lan Feb 26 '23 at 16:50
  • 1
    @Xiaoyan Lan can you manually type out (all or a subset of) the Ch_points and include it in the question? I can't seem to open the files in the link. – Ori Yarden PhD Feb 26 '23 at 17:24
  • @OriYarden I have tried to upload the subset of the dataset. However, I need to say that if there are just 50 records in Switzerland, it is easy to show all of them. After all, there are more than 15000 Swiss records in the original dataset. Meanwhile, around 10000 records are clustered in three detailed areas. I guess this is why it cannot work directly at the zoom level set by 'contextily'. – Xiaoyan Lan Feb 26 '23 at 20:01
  • The 'contextily.providers.OpenStreetMap.CH' has bound: [[45, 5], [48, 11]] which ranges from Lyon, France to Germany, and knowing that contextily can be finicky, perhaps providing a more specific source url could resolve this, such as 'url=image.getMapId(params)["tile_fetcher"].url_format' using import ee. – Ori Yarden PhD Feb 26 '23 at 20:52
  • @OriYarden Thanks a lot, Ori. I will use another specific URL to replace this one and try again. Have a nice day. – Xiaoyan Lan Feb 27 '23 at 15:21

0 Answers0