Hello everyone,
I am plotting my results through following code. I want to plot it with shapefile, with high quality map. Can you please respond how I add shapefile in this code and plot high quality map.
import pygmt
import pandas as pd
earthquake_data = pd.read_csv('catalog.csv')
stationfile="stations_BB.txt"
# Set up the figure
fig = pygmt.Figure()
# Define the region of interest
region = [89, 93, 26, 29]
# Set up the basemap with coastlines, rivers, and political boundaries
fig.basemap(region=region, projection="M6i", frame=True)
fig.coast(shorelines=True, resolution='10m', water='lightblue', land='lightgreen', borders=1)
# Plot earthquake data using circles (adjust size and color as needed)
fig.plot(
x=earthquake_data['longitude'],
y=earthquake_data['latitude'],
style='c0.1c', # Circle with a diameter of 0.1 cm
color='red',
label='Earthquakes'
)
# Add a legend
fig.legend(position='JMR+o0.2c', box=True)
# Add a title
fig.text(
text="Earthquake Catalog in Bhutan Region",
x=0.5,
y=1.05,
justify='CM',
font="12p,Helvetica-Bold",
)
# Save the figure to a file or display it
fig.show()
Welcome to the forum. Thanks for posting the PyGMT code that you are using. This example from the PyGMT Gallery can be a place to start with plotting data from shapefiles:
You need to use GeoPandas to read the shapefile format.
welcome to the GMT forum and sorry for my delayed response!
Currently I have difficulties to understand your issue: What do you mean by “plot it with shapefile, with high quality map”? Maybe it helps if you provide the shapefile you want to use.
As already mentioned by @EJFielding you can read a shapefile by using the Python package GeoPandas, for details please see the documentation at https://geopandas.org/en/stable/docs/reference/api/geopandas.read_file.html.
The data is read into a GeoDataFrame and the next steps (regarding plotting) depend on the kind of geometry (points, lines, or polygons). There are two PyGMT gallery examples which may serve as an orientation:
@yvonnefroehlich Thanks for your response. I have plot my catalog data with simple code. I want to add the country shape file with stations to plot my catalog inside the (https://data.humdata.org/dataset/cod-ab-btn?) the shapefile. I am facing an issue how I add this shape file in the code to plot my catalog inside that.
For completness: In your posted code you set resolution="10m". To get a higher resolution you may want to use resolution="h" (high):
import pygmt
# Set up the figure
fig = pygmt.Figure()
# Set up the basemap with coastlines, rivers, and political boundaries
fig.basemap(region=[89, 93, 26, 29], projection="M6i", frame=True)
fig.coast(shorelines=True, resolution="h", water="lightblue", land="lightgreen", borders=1)
# Save the figure to a file or display it
fig.show()
# fig.savefig(fname="plot_border_higherresolution.png")
@yvonnefroehlich Thank you very much for your detailed response. I have download the shapefile from the link (https://data.humdata.org/dataset/cod-ab-btn?) and placed it in the same foder where I have saved the code. I have run the code but it shows the bellow error. Please how I fix this error.
self.session.start(self, **kwargs)
File “fiona/ogrext.pyx”, line 588, in fiona.ogrext.Session.start
File “fiona/ogrext.pyx”, line 143, in fiona.ogrext.gdal_open_vector
fiona.errors.DriverError: btn_admbnda_bnlc_20201026_SHP/btn_admbnda_adm0_bnlc_20201026.shp: No such file or directory
Yes I have extract the zip file after downloading it, I have also changed the file path. The updated code which I am using is bellow but it shows again same error.
import pygmt
import geopandas as gpd
# Load shapefile into GeoDataFrame
# Different levels of the data are provided: adm0, adm1, adm2
# Please select for your needs
level = "adm0"
shapefile_path='/home/zamir/Downloads/gpd_btn/btn_admbnda_bnlc_20201026_SHP/btn_admbnda_adm0_bnlc_20201026.shp'
gpd_btn = gpd.read_file(shapefile_path)
gpd_btn = gpd.read_file(
f"{shapefile_path}btn_admbnda_bnlc_20201026_SHP/btn_admbnda_{level}_bnlc_20201026.shp"
)
# Set up new Figure instance
fig = pygmt.Figure()
fig.coast(
region=[89, 93, 26, 29],
projection="M10c",
land="lightgray",
frame=True,
)
# Plot the polygon(s) stored in the GeoDataFrame
fig.plot(
data=gpd_btn,
pen="0.5p,darkgreen",
fill="lightgreen",
)
fig.show()
# fig.savefig(fname=f"plot_shapefile_{level}.png")
The not found error is raised because the file path is built wrongly. In your code, you repeat the part "btn_admbnda_bnlc_20201026_SHP/btn_admbnda_adm0_bnlc_20201026.shp":