Repository: GenericMappingTools/pygmt

Issue:
Implement subplot command
This is a very cool new command.

Discussion:
willschlitzer: @weiji14 Any thought to using `sca` in a context manager? I know it doesn't have a `sca_end` equivalent the way that `inset` and `subplot` do, but I think using it with `with` statements will make it more clear that the code is creating different subplots, as opposed to a long line of similar looking plotting functions.
seisman: @weiji14 It's really a nice job! This is almost what I expect, although I haven't read the details.

A quick question about `ax` and `sca`. The long-form name of **-c** option is **panel** in GMT (see https://github.com/GenericMappingTools/gmt/blob/master/src/gmt_common_longoptions.h). Do we want to follow GMT or matplotlib?
weiji14: > @weiji14 Any thought to using `sca` in a context manager? I know it doesn't have a `sca_end` equivalent the way that `inset` and `subplot` do, but I think using it with `with` statements will make it more clear that the code is creating different subplots, as opposed to a long line of similar looking plotting functions.

I did think of this, but it's technically quite  hard to implement. What should a `with fig.sca()` block reset to? The first subplot? The next subplot? There might be a way using function decorators (?) to insert `ax=1` for all `fig.*` functions under the with block but I'll need more time to confirm that it's possible. Or feel free to suggest ideas :)

> @weiji14 It's really a nice job! This is almost what I expect, although I haven't read the details.

Thanks, it's really looking at @willschlitzer's `inset` wrapper that gave me the inspiration!

> A quick question about `ax` and `sca`. The long-form name of **-c** option is **panel** in GMT (see https://github.com/GenericMappingTools/gmt/blob/master/src/gmt_common_longoptions.h). Do we want to follow GMT or matplotlib?

Tough question. I think `panel` does seem like a name that conveys what is being done. But `ax` is shorter and more familiar to Python users :/ If only we could do both! Check with @GenericMappingTools/python and @willschlitzer.

seisman: Quick comments. In the 2nd example of the tutorial, you use:
```
fig.basemap(region=[0, 10, 0, 10], projection="X?", frame=["af", "WSne"], ax=axs[0, 0])
```

It seems passing the *index*, or *row,col* also works, 
```
fig.basemap(region=[0, 10, 0, 10], projection="X?", frame=["af", "WSne"], ax=0)
fig.basemap(region=[0, 10, 0, 10], projection="X?", frame=["af", "WSne"], ax="0,0")
```
but passing `[row, col]` doesn't:
```
fig.basemap(region=[0, 10, 0, 10], projection="X?", frame=["af", "WSne"], ax=[0, 0])
```
I think we can just set `ax` as the alias of `-c`, and convert `[0, 0]` to `0,0`, similar to what we're doing to `region` (converting `[0, 10, 0, 10]` to `0/10/0/10`.

BTW, it's still unclear to me how `ax=0` works in `basemap()`, as I don't see how you pass **-c** to `basemap`, or calling `subplot set`. Could you please explain it a little bit?

weiji14: > Quick comments. In the 2nd example of the tutorial, you use:
> 
> ```
> fig.basemap(region=[0, 10, 0, 10], projection="X?", frame=["af", "WSne"], ax=axs[0, 0])
> ```
> 
> It seems passing the _index_, or _row,col_ also works,
> 
> ```
> fig.basemap(region=[0, 10, 0, 10], projection="X?", frame=["af", "WSne"], ax=0)
> fig.basemap(region=[0, 10, 0, 10], projection="X?", frame=["af", "WSne"], ax="0,0")
> ```

Yes, in pure GMT, you would use something like `-c0` or `-c0,1`. I'm debating whether to even have the `axs` numpy array thing, it's meant to be like matplotlib's `Axes` object but we can't do `ax.plot` so maybe that will confuse users?

> but passing `[row, col]` doesn't:
> 
> ```
> fig.basemap(region=[0, 10, 0, 10], projection="X?", frame=["af", "WSne"], ax=[0, 0])
> ```
> 
> I think we can just set `ax` as the alias of `-c`, and convert `[0, 0]` to `0,0`, similar to what we're doing to `region` (converting `[0, 10, 0, 10]` to `0/10/0/10`.

Ok, I'll read this as aliasing `c` to `ax` (rather than using `panel`), and we just need a `sequence (slash)` converter in place.

> BTW, it's still unclear to me how `ax=0` works in `basemap()`, as I don't see how you pass **-c** to `basemap`, or calling `subplot set`. Could you please explain it a little bit?

`ax` is an alias for `c`, so `ax=0` is equivalent to `-c0` in pure GMT. Using `fig.basemap(..., ax=0)` would tell the program to directly plot in the first (0) subplot, so `subplot set`/`fig.sca()` is not needed. See also https://docs.generic-mapping-tools.org/latest/gmt.html#c-full.

seisman: > `ax` is an alias for `c`, so `ax=0` is equivalent to `-c0` in pure GMT. Using `fig.basemap(..., ax=0)` would tell the program to directly plot in the first (0) subplot, so `subplot set`/`fig.sca()` is not needed. See also https://docs.generic-mapping-tools.org/latest/gmt.html#c-full.

I see. I was fooled by the github web UI. When I review the changes, I thought you make the alias (`ax`<=>`c`) in `contour()` but can pass `ax` to `basemap()`. What magic! Now it all make sense to me.

> Yes, in pure GMT, you would use something like `-c0` or `-c0,1`. I'm debating whether to even have the `axs` numpy array thing, it's meant to be like matplotlib's `Axes` object but we can't do `ax.plot` so maybe that will confuse users?

Yes, especially that if `ax = [0, 0]` works, why do we bother to type `ax=axs[0, 0]`.

willschlitzer: > 
> Yeah, but map elements outside a with fig.sca() block can't just go into the main figure, as it would need to be outside the       
> fig.subplot() indent. It's hard to come up with a good solution here.

I hope I didn't cause too much trouble with suggesting a context manager!

> So should we 1) just go back to having fig.sca as an ordinary function (i.e. not a context manager)? Or 2) keep fig.sca() as a 
> context manager so that people can organize their subplots using indentation (if they choose to, or not)?

I'm a fan of option 2. While I understand that `fig.sca()` doesn't need to be "closed" the way that `inset` and `subplot` do, I still think it makes it visually easier to see different `with` statements being used to plot on different subplots. The good news here is that if someone forgets to indent a line underneath a context manager it will still plot to the appropriate subplot figure, so it's not too punishing if they have a syntax typo.

weiji14: > > Yeah, but map elements outside a with fig.sca() block can't just go into the main figure, as it would need to be outside the
> > fig.subplot() indent. It's hard to come up with a good solution here.
> 
> I hope I didn't cause too much trouble with suggesting a context manager!

It's actually a good idea, just that we're technically limited by GMT! We're just trying to find a way that's less confusing than matplotlib if possible.

> > So should we 1) just go back to having fig.sca as an ordinary function (i.e. not a context manager)? Or 2) keep fig.sca() as a
> > context manager so that people can organize their subplots using indentation (if they choose to, or not)?
> 
> I'm a fan of option 2. While I understand that `fig.sca()` doesn't need to be "closed" the way that `inset` and `subplot` do, I still think it makes it visually easier to see different `with` statements being used to plot on different subplots. The good news here is that if someone forgets to indent a line underneath a context manager it will still plot to the appropriate subplot figure, so it's not too punishing if they have a syntax typo.

I'm inclined to agree with @willschlitzer here, but after re-reading through the [Zen of Python](https://www.python.org/dev/peps/pep-0020) (my favourite is "There should be one-- and preferably only one --obvious way to do it."), I'm a bit split on what to choose. Having the `with` block option is nice as it allows for more flexibility, but is it too much freedom?

One thing I'll point out too is that the 'many' ways of making a subplot depends on how fancy people want to go, and maybe we could just document it better in the tutorial:

0. Use `fig.shift_origin()` #289.

1. The `ax` way for fairly simple subplots, E.g.:

```python
fig.basemap(..., ax=[0, 1])  # This is in subplot 1
fig.coast(...)  # This is also placed in subplot 1
fig.basemap(..., ax=[0, 1])  # This is placed in subplot 2
```

2. The `with fig.sca(ax=...):` way, better for complicated subplots if you want to override the title. E.g.:

```python
with fig.sca(ax=0, fixedlabel="subplot 0"):
    fig.basemap(...)
with fig.sca(ax=1, fixedlabel="subplot 1"):
    fig.basemap(...)
```
willschlitzer: @weiji14 Not to throw a last-minute change right before we're also trying to get this PR closed out, but did you ever look at making a Subplot() class that is a sub-class of Figure(), with every subplot being it's own Python object? I realize this is a complete change over what I suggested previously, but I think that would provide a clear visual breakdown of which plotting elements are applied to which subplot without needing a context manager/with statement.
weiji14: > @weiji14 Not to throw a last-minute change right before we're also trying to get this PR closed out, but did you ever look at making a Subplot() class that is a sub-class of Figure(), with every subplot being it's own Python object?

Yes, this (Subplot() sub-class) was pretty much the previous v2 design at #427, though not exactly if I follow your reasoning. One problem with that previous design was that everything in a Subplot instance was a subplot, and you couldn't add say, a colorbar or legend outside of the subplot for the whole figure (after `subplot end`). The thing is that a 'Subplot' (all panels) exists **inside** of a larger 'Figure'.

> I realize this is a complete change over what I suggested previously, but I think that would provide a clear visual breakdown of which plotting elements are applied to which subplot without needing a context manager/with statement.

Actually, could you make a rough skeleton code sketch of how this would work? I'm actually not good with Python classes, more of a functional programming type :blush:
willschlitzer: I'll take a look at it! I'm not too smart about Python classes as well, so I could be thinking of class inheritance that doesn't actually work.

For the sake of getting this wrapped, since I would consider it the flagship new function of v0.3, I like the function as it currently stands and don't think it's bad to have multiple ways to create subplots, as long as we create sufficient tutorials for it.
weiji14: > For the sake of getting this wrapped, since I would consider it the flagship new function of v0.3, I like the function as it currently stands and don't think it's bad to have multiple ways to create subplots, as long as we create sufficient tutorials for it.

Yes, it's been a long time coming since the original implementation in #412 (Apr 2020), and I really want to close our oldest issue #20 (Jul 2017)! 
willschlitzer: I'm waiting for the tests to finish running and I'll squash and merge #808 in the next few minutes.
weiji14: Any other comments @GenericMappingTools/python-maintainers? I'd like to get this done by the weekend if possible (but don't hesitate to ping me on any quality issues). Check the following documentation pages and see if there's any mistakes:

- https://pygmt-git-subplotv3.gmt.vercel.app/tutorials/subplots.html
- https://pygmt-git-subplotv3.gmt.vercel.app/api/generated/pygmt.Figure.subplot.html

Patch:
diff --git a/doc/api/index.rst b/doc/api/index.rst
index e6b45603..0ef0b359 100644
--- a/doc/api/index.rst
+++ b/doc/api/index.rst
@@ -36,7 +36,9 @@ Plotting data and laying out the map:
     Figure.meca
     Figure.plot
     Figure.plot3d
+    Figure.set_panel
     Figure.shift_origin
+    Figure.subplot
     Figure.text
 
 Color palette table generation:
diff --git a/doc/index.rst b/doc/index.rst
index 399c7ad1..7dcc11b8 100644
--- a/doc/index.rst
+++ b/doc/index.rst
@@ -39,6 +39,7 @@
     tutorials/earth-relief.rst
     tutorials/3d-perspective-image.rst
     tutorials/inset.rst
+    tutorials/subplots.rst
     tutorials/configuration.rst
 
 .. toctree::
diff --git a/examples/tutorials/subplots.py b/examples/tutorials/subplots.py
new file mode 100644
index 00000000..64c5e045
--- /dev/null
+++ b/examples/tutorials/subplots.py
@@ -0,0 +1,237 @@
+"""
+Making subplots
+===============
+
+When you're preparing a figure for a paper, there will often be times when
+you'll need to put many individual plots into one large figure, and label them
+'abcd'. These individual plots are called subplots.
+
+There are two main ways to create subplots in GMT:
+
+- Use :meth:`pygmt.Figure.shift_origin` to manually move each individual plot
+  to the right position.
+- Use :meth:`pygmt.Figure.subplot` to define the layout of the subplots.
+
+The first method is easier to use and should handle simple cases involving a
+couple of subplots. For more advanced subplot layouts, however, we recommend the
+use of :meth:`pygmt.Figure.subplot` which offers finer grained control, and
+this is what the tutorial below will cover.
+"""
+# sphinx_gallery_thumbnail_number = 3
+
+import pygmt
+
+###############################################################################
+#
+# Let's start by initializing a :class:`pygmt.Figure` instance.
+
+fig = pygmt.Figure()
+
+###############################################################################
+# Define subplot layout
+# ---------------------
+#
+# The :meth:`pygmt.Figure.subplot` function is used to set up the layout, size,
+# and other attributes of the figure. It divides the whole canvas into regular
+# grid areas with *n* rows and *m* columns. Each grid area can contain an
+# individual subplot. For example:
+
+###############################################################################
+# .. code-block:: default
+#
+#     with fig.subplot(nrows=2, ncols=3, figsize=("15c", "6c"), frame="lrtb"):
+#         ...
+
+###############################################################################
+# will define our figure to have a 2 row and 3 column grid layout.
+# ``figsize=("15c", "6c")`` defines the overall size of the figure to be 15 cm
+# wide by 6 cm high. Using ``frame="lrtb"`` allows us to customize the map frame
+# for all subplots instead of setting them individually. The figure layout will
+# look like the following:
+
+with fig.subplot(nrows=2, ncols=3, figsize=("15c", "6c"), frame="lrtb"):
+    for i in range(2):  # row number starting from 0
+        for j in range(3):  # column number starting from 0
+            index = i * 3 + j  # index number starting from 0
+            with fig.set_panel(panel=index):  # sets the current panel
+                fig.text(
+                    position="MC",
+                    text=f"index: {index}; row: {i}, col: {j}",
+                    region=[0, 1, 0, 1],
+                )
+fig.show()
+
+###############################################################################
+# The :meth:`pygmt.Figure.set_panel` function activates a specified subplot, and
+# all subsequent plotting functions will take place in that subplot panel. This
+# is similar to matplotlib's ``plt.sca`` method. In order to specify a subplot,
+# you will need to provide the identifier for that subplot via the ``panel``
+# parameter. Pass in either the *index* number, or a tuple/list like
+# (*row*, *col*) to ``panel``.
+
+###############################################################################
+# .. note::
+#
+#     The row and column numbering starts from 0. So for a subplot layout with
+#     N rows and M columns, row numbers will go from 0 to N-1, and column
+#     numbers will go from 0 to M-1.
+
+###############################################################################
+# For example, to activate the subplot on the top right corner (index: 2) at
+# *row*\=0 and *col*\=2, so that all subsequent plotting commands happen
+# there, you can use the following command:
+
+###############################################################################
+# .. code-block:: default
+#
+#     with fig.set_panel(panel=[0, 2]):
+#         ...
+
+###############################################################################
+# Making your first subplot
+# -------------------------
+# Next, let's use what we learned above to make a 2 row by 2 column subplot
+# figure. We'll also pick up on some new parameters to configure our subplot.
+
+fig = pygmt.Figure()
+with fig.subplot(
+    nrows=2,
+    ncols=2,
+    figsize=("15c", "6c"),
+    autolabel=True,
+    frame=["af", "WSne"],
+    margins=["0.1c", "0.2c"],
+    title="My Subplot Heading",
+):
+    fig.basemap(region=[0, 10, 0, 10], projection="X?", panel=[0, 0])
+    fig.basemap(region=[0, 20, 0, 10], projection="X?", panel=[0, 1])
+    fig.basemap(region=[0, 10, 0, 20], projection="X?", panel=[1, 0])
+    fig.basemap(region=[0, 20, 0, 20], projection="X?", panel=[1, 1])
+fig.show()
+
+###############################################################################
+# In this example, we define a 2-row, 2-column (2x2) subplot layout using
+# :meth:`pygmt.Figure.subplot`. The overall figure dimensions is set to be
+# 15 cm wide and 6 cm high (``figsize=["15c", "6c"]``). In addition, we use
+# some optional parameters to fine-tune some details of the figure creation:
+#
+# - ``autolabel=True``: Each subplot is automatically labelled abcd
+# - ``margins=["0.1c", "0.2c"]``: adjusts the space between adjacent subplots.
+#   In this case, it is set as 0.1 cm in the X direction and 0.2 cm in the Y
+#   direction.
+# - ``title="My Subplot Heading"``: adds a title on top of the whole figure.
+#
+# Notice that each subplot was set to use a linear projection ``"X?"``.
+# Usually, we need to specify the width and height of the map frame, but it is
+# also possible to use a question mark ``"?"`` to let GMT decide automatically
+# on what is the most appropriate width/height for the each subplot's map
+# frame.
+
+###############################################################################
+# .. tip::
+#
+#     In the above example, we used the following commands to activate the
+#     four subplots explicitly one after another::
+#
+#         fig.basemap(..., panel=[0, 0])
+#         fig.basemap(..., panel=[0, 1])
+#         fig.basemap(..., panel=[1, 0])
+#         fig.basemap(..., panel=[1, 1])
+#
+#     In fact, we can just use ``fig.basemap(..., panel=True)`` without
+#     specifying any subplot index number, and GMT will automatically activate
+#     the next subplot panel.
+
+###############################################################################
+# .. note::
+#
+#     All plotting functions (e.g. :meth:`pygmt.Figure.coast`,
+#     :meth:`pygmt.Figure.text`, etc) are able to use ``panel`` parameter when
+#     in subplot mode. Once a panel is activated using ``panel`` or
+#     :meth:`pygmt.Figure.set_panel`, subsequent plotting commands that don't
+#     set a ``panel`` will have their elements added to the same panel as
+#     before.
+
+###############################################################################
+# Shared X and Y axis labels
+# --------------------------
+# In the example above with the four subplots, the two subplots for each row
+# have the same Y-axis range, and the two subplots for each column have the
+# same X-axis range. You can use the ``sharex``/``sharey`` parameters to set a
+# common X and/or Y axis between subplots.
+
+fig = pygmt.Figure()
+with fig.subplot(
+    nrows=2,
+    ncols=2,
+    figsize=("15c", "6c"),  # width of 15 cm, height of 6 cm
+    autolabel=True,
+    margins=["0.3c", "0.2c"],  # horizontal 0.3 cm and vertical 0.2 cm margins
+    title="My Subplot Heading",
+    sharex="b",  # shared x-axis on the bottom side
+    sharey="l",  # shared y-axis on the left side
+    frame="WSrt",
+):
+    fig.basemap(region=[0, 10, 0, 10], projection="X?", panel=True)
+    fig.basemap(region=[0, 20, 0, 10], projection="X?", panel=True)
+    fig.basemap(region=[0, 10, 0, 20], projection="X?", panel=True)
+    fig.basemap(region=[0, 20, 0, 20], projection="X?", panel=True)
+fig.show()
+
+###############################################################################
+# ``sharex="b"`` indicates that subplots in a column will share the x-axis, and
+# only the **b**\ ottom axis is displayed. ``sharey="l"`` indicates that
+# subplots within a row will share the y-axis, and only the **l**\ eft axis is
+# displayed.
+#
+# Of course, instead of using the ``sharex``/``sharey`` option, you can also
+# set a different ``frame`` for each subplot to control the axis properties
+# individually for each subplot.
+
+###############################################################################
+# Advanced subplot layouts
+# ------------------------
+#
+# Nested subplot are currently not supported. If you want to create more
+# complex subplot layouts, some manual adjustments are needed.
+#
+# The following example draws three subplots in a 2-row, 2-column layout, with
+# the first subplot occupying the first row.
+
+fig = pygmt.Figure()
+# Bottom row, two subplots
+with fig.subplot(nrows=1, ncols=2, figsize=("15c", "3c"), autolabel="b)"):
+    fig.basemap(
+        region=[0, 5, 0, 5], projection="X?", frame=["af", "WSne"], panel=[0, 0]
+    )
+    fig.basemap(
+        region=[0, 5, 0, 5], projection="X?", frame=["af", "WSne"], panel=[0, 1]
+    )
+# Move plot origin by 1 cm above the height of the entire figure
+fig.shift_origin(yshift="h+1c")
+# Top row, one subplot
+with fig.subplot(nrows=1, ncols=1, figsize=("15c", "3c"), autolabel="a)"):
+    fig.basemap(
+        region=[0, 10, 0, 10], projection="X?", frame=["af", "WSne"], panel=[0, 0]
+    )
+    fig.text(text="TEXT", x=5, y=5)
+
+fig.show()
+
+###############################################################################
+#
+# We start by drawing the bottom two subplots, setting ``autolabel="b)"`` so
+# that the subplots are labelled 'b)' and 'c)'. Next, we use
+# :meth:`pygmt.Figure.shift_origin` to move the plot origin 1 cm above the
+# **h**\ eight of the entire figure that is currently plotted (i.e. the bottom
+# row subplots). A single subplot is then plotted on the top row. You may need
+# to adjust the ``yshift`` parameter to make your plot look nice. This top row
+# uses ``autolabel="a)"``, and we also plotted some text inside. Note that
+# ``projection="X?"`` was used to let GMT automatically determine the size of
+# the subplot according to the size of the subplot area.
+
+###############################################################################
+# You can also manually override the ``autolabel`` for each subplot using for
+# example, ``fig.set_panel(..., fixedlabel="b) Panel 2")`` which would allow
+# you to manually label a single subplot as you wish. This can be useful for
+# adding a more descriptive subtitle to individual subplots.
diff --git a/pygmt/figure.py b/pygmt/figure.py
index 09ce2223..ba43963f 100644
--- a/pygmt/figure.py
+++ b/pygmt/figure.py
@@ -387,5 +387,7 @@ class Figure:
         meca,
         plot,
         plot3d,
+        set_panel,
+        subplot,
         text,
     )
diff --git a/pygmt/src/__init__.py b/pygmt/src/__init__.py
index 016e610d..3ef9af36 100644
--- a/pygmt/src/__init__.py
+++ b/pygmt/src/__init__.py
@@ -24,6 +24,7 @@ from pygmt.src.makecpt import makecpt
 from pygmt.src.meca import meca
 from pygmt.src.plot import plot
 from pygmt.src.plot3d import plot3d
+from pygmt.src.subplot import set_panel, subplot
 from pygmt.src.surface import surface
 from pygmt.src.text import text_ as text  # "text" is an argument within "text_"
 from pygmt.src.which import which
diff --git a/pygmt/src/subplot.py b/pygmt/src/subplot.py
new file mode 100644
index 00000000..fdf11856
--- /dev/null
+++ b/pygmt/src/subplot.py
@@ -0,0 +1,232 @@
+"""
+subplot - Manage modern mode figure subplot configuration and selection.
+"""
+import contextlib
+
+from pygmt.clib import Session
+from pygmt.exceptions import GMTInvalidInput
+from pygmt.helpers import (
+    build_arg_string,
+    fmt_docstring,
+    is_nonstr_iter,
+    kwargs_to_strings,
+    use_alias,
+)
+
+
+@fmt_docstring
+@contextlib.contextmanager
+@use_alias(
+    Ff="figsize",
+    Fs="subsize",
+    A="autolabel",
+    B="frame",
+    C="clearance",
+    J="projecton",
+    M="margins",
+    R="region",
+    SC="sharex",
+    SR="sharey",
+    T="title",
+    V="verbose",
+    X="xshift",
+    Y="yshift",
+)
+@kwargs_to_strings(Ff="sequence", Fs="sequence", M="sequence", R="sequence")
+def subplot(self, nrows=1, ncols=1, **kwargs):
+    r"""
+    Create multi-panel subplot figures.
+
+    This function is used to split the current figure into a rectangular layout
+    of subplots that each may contain a single self-contained figure. Begin by
+    defining the layout of the entire multi-panel illustration. Several
+    parameters are available to specify the systematic layout, labeling,
+    dimensions, and more for the subplots.
+
+    Full option list at :gmt-docs:`subplot.html#synopsis-begin-mode`
+
+    {aliases}
+
+    Parameters
+    ----------
+    nrows : int
+        Number of vertical rows of the subplot grid.
+    ncols : int
+        Number of horizontal columns of the subplot grid.
+    figsize : tuple
+        Specify the final figure dimensions as (*width*, *height*).
+    subsize : tuple
+        Specify the dimensions of each subplot directly as (*width*, *height*).
+        Note that only one of ``figsize`` or ``subsize`` can be provided at
+        once.
+
+    autolabel : bool or str
+        [*autolabel*][**+c**\ *dx*\ [/*dy*]][**+g**\ *fill*][**+j**\|\ **J**\
+        *refpoint*][**+o**\ *dx*\ [/*dy*]][**+p**\ *pen*][**+r**\|\ **R**]
+        [**+v**].
+        Specify automatic tagging of each subplot. Append either a number or
+        letter [a]. This sets the tag of the first, top-left subplot and others
+        follow sequentially. Surround the number or letter by parentheses on
+        any side if these should be typeset as part of the tag. Use
+        **+j**\|\ **J**\ *refpoint* to specify where the tag should be placed
+        in the subplot [TL]. Note: **+j** sets the justification of the tag to
+        *refpoint* (suitable for interior tags) while **+J** instead selects
+        the mirror opposite (suitable for exterior tags). Append
+        **+c**\ *dx*\[/*dy*] to set the clearance between the tag and a
+        surrounding text box requested via **+g** or **+p** [3p/3p, i.e., 15%
+        of the :gmt-term:`FONT_TAG` size dimension]. Append **+g**\ *fill* to
+        paint the tag's text box with *fill* [no painting]. Append
+        **+o**\ *dx*\ [/*dy*] to offset the tag's reference point in the
+        direction implied by the justification [4p/4p, i.e., 20% of the
+        :gmt-term:`FONT_TAG` size]. Append **+p**\ *pen* to draw the outline of
+        the tag's text box using selected *pen* [no outline]. Append **+r** to
+        typeset your tag numbers using lowercase Roman numerals; use **+R** for
+        uppercase Roman numerals [Arabic numerals]. Append **+v** to increase
+        tag numbers vertically down columns [horizontally across rows].
+    {B}
+    clearance : str or list
+        [*side*]\ *clearance*.
+        Reserve a space of dimension *clearance* between the margin and the
+        subplot on the specified side, using *side* values from **w**, **e**,
+        **s**, or **n**; or **x** for both **w** and **e**; or **y** for both
+        **s** and **n**. No *side* means all sides (i.e. ``clearance='1c'``
+        would set a clearance of 1 cm on all sides). The option is repeatable
+        to set aside space on more than one side (e.g. ``clearance=['w1c',
+        's2c']`` would set a clearance of 1 cm on west side and 2 cm on south
+        side). Such space will be left untouched by the main map plotting but
+        can be accessed by modules that plot scales, bars, text, etc.
+    {J}
+    margins : str or list
+        This is margin space that is added between neighboring subplots (i.e.,
+        the interior margins) in addition to the automatic space added for tick
+        marks, annotations, and labels. The margins can be specified as either:
+
+        - a single value (for same margin on all sides). E.g. '5c'.
+        - a pair of values (for setting separate horizontal and vertical
+          margins). E.g. ['5c', '3c'].
+        - a set of four values (for setting separate left, right, bottom, and
+          top margins). E.g. ['1c', '2c', '3c', '4c'].
+
+        The actual gap created is always a sum of the margins for the two
+        opposing sides (e.g., east plus west or south plus north margins)
+        [Default is half the primary annotation font size, giving the full
+        annotation font size as the default gap].
+    {R}
+    sharex : bool or str
+        Set subplot layout for shared x-axes. Use when all subplots in a column
+        share a common *x*-range. If ``sharex=True``, the first (i.e.,
+        **t**\ op) and the last (i.e., **b**\ ottom) rows will have
+        *x*-annotations; use ``sharex='t'`` or ``sharex='b'`` to select only
+        one of those two rows [both]. Append **+l** if annotated *x*-axes
+        should have a label [none]; optionally append the label if it is the
+        same for the entire subplot. Append **+t** to make space for subplot
+        titles for each row; use **+tc** for top row titles only [no subplot
+        titles].
+    sharey : bool or str
+        Set subplot layout for shared y-axes. Use when all subplots in a row
+        share a common *y*-range. If ``sharey=True``, the first (i.e.,
+        **l**\ eft) and the last (i.e., **r**\ ight) columns will have
+        *y*-annotations; use ``sharey='l'`` or ``sharey='r'`` to select only
+        one of those two columns [both]. Append **+l** if annotated *y*-axes
+        will have a label [none]; optionally, append the label if it is the
+        same for the entire subplot. Append **+p** to make all annotations
+        axis-parallel [horizontal]; if not used you may have to set
+        ``clearance`` to secure extra space for long horizontal annotations.
+
+        Notes for ``sharex``/``sharey``:
+
+        - Labels and titles that depends on which row or column are specified
+          as usual via a subplot's own ``frame`` setting.
+        - Append **+w** to the ``figsize`` or ``subsize`` parameter to draw
+          horizontal and vertical lines between interior panels using selected
+          pen [no lines].
+    title : str
+        While individual subplots can have titles (see ``sharex``/``sharey`` or
+        ``frame``), the entire figure may also have an overarching *heading*
+        [no heading]. Font is determined by setting :gmt-term:`FONT_HEADING`.
+    {V}
+    {XY}
+    """
+    kwargs = self._preprocess(**kwargs)  # pylint: disable=protected-access
+    # allow for spaces in string with needing double quotes
+    kwargs["A"] = f'"{kwargs.get("A")}"' if kwargs.get("A") is not None else None
+    kwargs["T"] = f'"{kwargs.get("T")}"' if kwargs.get("T") else None
+
+    if nrows < 1 or ncols < 1:
+        raise GMTInvalidInput("Please ensure that both 'nrows'>=1 and 'ncols'>=1.")
+    if kwargs.get("Ff") and kwargs.get("Fs"):
+        raise GMTInvalidInput(
+            "Please provide either one of 'figsize' or 'subsize' only."
+        )
+
+    with Session() as lib:
+        try:
+            arg_str = " ".join(["begin", f"{nrows}x{ncols}", build_arg_string(kwargs)])
+            lib.call_module("subplot", arg_str)
+            yield
+        finally:
+            v_arg = build_arg_string({"V": kwargs.get("V")})
+            lib.call_module("subplot", f"end {v_arg}")
+
+
+@fmt_docstring
+@contextlib.contextmanager
+@use_alias(A="fixedlabel", C="clearance", V="verbose")
+def set_panel(self, panel=None, **kwargs):
+    r"""
+    Set the current subplot panel to plot on.
+
+    Before you start plotting you must first select the active subplot. Note:
+    If any *projection* option is passed with the question mark **?** as scale
+    or width when plotting subplots, then the dimensions of the map are
+    automatically determined by the subplot size and your region. For Cartesian
+    plots: If you want the scale to apply equally to both dimensions then you
+    must specify ``projection="x"`` [The default ``projection="X"`` will fill
+    the subplot by using unequal scales].
+
+    {aliases}
+
+    Parameters
+    ----------
+    panel : str or list
+        *row,col*\|\ *index*.
+        Sets the current subplot until further notice. **Note**: First *row*
+        or *col* is 0, not 1. If not given we go to the next subplot by order
+        specified via ``autolabel`` in :meth:`pygmt.Figure.subplot`. As an
+        alternative, you may bypass using :meth:`pygmt.Figure.set_panel` and
+        instead supply the common option **panel**\ =[*row,col*] to the first
+        plot command you issue in that subplot. GMT maintains information about
+        the current figure and subplot. Also, you may give the one-dimensional
+        *index* instead which starts at 0 and follows the row or column order
+        set via ``autolabel`` in :meth:`pygmt.Figure.subplot`.
+
+    fixedlabel : str
+        Overrides the automatic labeling with the given string. No modifiers
+        are allowed. Placement, justification, etc. are all inherited from how
+        ``autolabel`` was specified by the initial :meth:`pygmt.Figure.subplot`
+        command.
+
+    clearance : str or list
+        [*side*]\ *clearance*.
+        Reserve a space of dimension *clearance* between the margin and the
+        subplot on the specified side, using *side* values from **w**, **e**,
+        **s**, or **n**. The option is repeatable to set aside space on more
+        than one side (e.g. ``clearance=['w1c', 's2c']`` would set a clearance
+        of 1 cm on west side and 2 cm on south side). Such space will be left
+        untouched by the main map plotting but can be accessed by modules that
+        plot scales, bars, text, etc. This setting overrides the common
+        clearances set by ``clearance`` in the initial
+        :meth:`pygmt.Figure.subplot` call.
+
+    {V}
+    """
+    kwargs = self._preprocess(**kwargs)  # pylint: disable=protected-access
+    # allow for spaces in string with needing double quotes
+    kwargs["A"] = f'"{kwargs.get("A")}"' if kwargs.get("A") is not None else None
+    # convert tuple or list to comma-separated str
+    panel = ",".join(map(str, panel)) if is_nonstr_iter(panel) else panel
+
+    with Session() as lib:
+        arg_str = " ".join(["set", f"{panel}", build_arg_string(kwargs)])
+        lib.call_module(module="subplot", args=arg_str)
+        yield


Test patch:
diff --git a/pygmt/tests/test_subplot.py b/pygmt/tests/test_subplot.py
new file mode 100644
index 00000000..44d78a6a
--- /dev/null
+++ b/pygmt/tests/test_subplot.py
@@ -0,0 +1,110 @@
+"""
+Tests subplot.
+"""
+import pytest
+from pygmt import Figure
+from pygmt.exceptions import GMTInvalidInput
+from pygmt.helpers.testing import check_figures_equal
+
+
+@check_figures_equal()
+def test_subplot_basic_frame():
+    """
+    Create a subplot figure with 1 vertical row and 2 horizontal columns, and
+    ensure map frame setting is applied to all subplot figures.
+    """
+    fig_ref, fig_test = Figure(), Figure()
+    with fig_ref.subplot(nrows=1, ncols=2, Ff="6c/3c", B="WSne"):
+        with fig_ref.set_panel(panel=0):
+            fig_ref.basemap(region=[0, 3, 0, 3], frame="+tplot0")
+        with fig_ref.set_panel(panel=1):
+            fig_ref.basemap(region=[0, 3, 0, 3], frame="+tplot1")
+    with fig_test.subplot(nrows=1, ncols=2, figsize=("6c", "3c"), frame="WSne"):
+        with fig_test.set_panel(panel="0,0"):
+            fig_test.basemap(region=[0, 3, 0, 3], frame="+tplot0")
+        with fig_test.set_panel(panel=[0, 1]):
+            fig_test.basemap(region=[0, 3, 0, 3], frame="+tplot1")
+    return fig_ref, fig_test
+
+
+@check_figures_equal()
+def test_subplot_direct():
+    """
+    Plot map elements to subplot directly using the panel parameter.
+    """
+    fig_ref, fig_test = Figure(), Figure()
+    with fig_ref.subplot(nrows=2, ncols=1, Fs="3c/3c"):
+        fig_ref.basemap(region=[0, 3, 0, 3], frame="af", panel=0)
+        fig_ref.basemap(region=[0, 3, 0, 3], frame="af", panel=1)
+    with fig_test.subplot(nrows=2, ncols=1, subsize=("3c", "3c")):
+        fig_test.basemap(region=[0, 3, 0, 3], frame="af", panel=[0, 0])
+        fig_test.basemap(region=[0, 3, 0, 3], frame="af", panel=[1, 0])
+    return fig_ref, fig_test
+
+
+@check_figures_equal()
+def test_subplot_autolabel_margins_title():
+    """
+    Make subplot figure with autolabels, setting some margins and a title.
+    """
+    fig_ref, fig_test = Figure(), Figure()
+    kwargs = dict(nrows=2, ncols=1, figsize=("15c", "6c"))
+
+    with fig_ref.subplot(A="a)", M="0.3c/0.1c", T="Subplot Title", **kwargs):
+        fig_ref.basemap(region=[0, 1, 2, 3], frame="WSne", c="0,0")
+        fig_ref.basemap(region=[4, 5, 6, 7], frame="WSne", c="1,0")
+
+    with fig_test.subplot(
+        autolabel=True, margins=["0.3c", "0.1c"], title="Subplot Title", **kwargs
+    ):
+        fig_test.basemap(region=[0, 1, 2, 3], frame="WSne", panel=[0, 0])
+        fig_test.basemap(region=[4, 5, 6, 7], frame="WSne", panel=[1, 0])
+
+    return fig_ref, fig_test
+
+
+@check_figures_equal()
+def test_subplot_clearance_and_shared_xy_axis_layout():
+    """
+    Ensure subplot clearance works, and that the layout can be set to use
+    shared X and Y axis labels across columns and rows.
+    """
+    fig_ref, fig_test = Figure(), Figure()
+    kwargs = dict(nrows=2, ncols=2, frame="WSrt", figsize=("5c", "5c"))
+
+    with fig_ref.subplot(C="y0.2c", SC="t", SR="", **kwargs):
+        fig_ref.basemap(region=[0, 4, 0, 4], projection="X?", panel=True)
+        fig_ref.basemap(region=[0, 8, 0, 4], projection="X?", panel=True)
+        fig_ref.basemap(region=[0, 4, 0, 8], projection="X?", panel=True)
+        fig_ref.basemap(region=[0, 8, 0, 8], projection="X?", panel=True)
+
+    with fig_test.subplot(
+        clearance=["s0.2c", "n0.2c"], sharex="t", sharey=True, **kwargs
+    ):
+        fig_test.basemap(region=[0, 4, 0, 4], projection="X?", panel=True)
+        fig_test.basemap(region=[0, 8, 0, 4], projection="X?", panel=True)
+        fig_test.basemap(region=[0, 4, 0, 8], projection="X?", panel=True)
+        fig_test.basemap(region=[0, 8, 0, 8], projection="X?", panel=True)
+
+    return fig_ref, fig_test
+
+
+def test_subplot_figsize_and_subsize_error():
+    """
+    Check that an error is raised when both figsize and subsize parameters are
+    passed into subplot.
+    """
+    fig = Figure()
+    with pytest.raises(GMTInvalidInput):
+        with fig.subplot(figsize=("2c", "1c"), subsize=("2c", "1c")):
+            pass
+
+
+def test_subplot_nrows_ncols_less_than_one_error():
+    """
+    Check that an error is raised when nrows or ncols is less than one.
+    """
+    fig = Figure()
+    with pytest.raises(GMTInvalidInput):
+        with fig.subplot(nrows=0, ncols=-1, figsize=("2c", "1c")):
+            pass