modflowpy / flopy

A Python package to create, run, and post-process MODFLOW-based models.
https://flopy.readthedocs.io
Other
517 stars 313 forks source link

PlotMapView resets limits of existing axes to entire model domain #2101

Closed rubencalje closed 8 months ago

rubencalje commented 8 months ago

Describe the bug Maybe this is intentional and not a bug, but: PlotMapView resets the matplotlib axes limits to the entire model domain.

To Reproduce Steps to reproduce the behavior: Run the following code:

import numpy as np
import flopy
import matplotlib.pyplot as plt

nrow=30
ncol=20
delr=np.array([100] * ncol)
delc = np.array([100] * nrow)
modelgrid = flopy.discretization.StructuredGrid(nrow=nrow, ncol=ncol, delr=delr, delc=delc)

f, ax = plt.subplots()
ax.axis('scaled')
ax.axis([500, 800, 200, 400])
pmv = flopy.plot.PlotMapView(modelgrid=modelgrid)
pmv.plot_grid()

Expected behavior I expect the ax to have an extent specified in ax.axis(): [500, 800, 200, 400]. However, the axes is zoomed out to fit the entire model domain.

Screenshots The outcome of the code above is: outcome

While the desired outcome would be: expected_outcome

Possible solution If you agree this is a bug, I can add a pull request, that checks if autoscaling is on when an axes is supplied as a parameter to PlotMapView, see https://matplotlib.org/stable/users/explain/axes/autoscale.html#controlling-autoscale

If autoscaling is on, the axes limits are still set to the entire model domain. However when autoscaling is off, the axes limits are not changed.

langevin-usgs commented 8 months ago

Hey Ruben, I think the easy way to handle this is to make the changes to ax later on. Maybe something like this?

f, ax = plt.subplots()
pmv = flopy.plot.PlotMapView(modelgrid=modelgrid, ax=ax)
pmv.plot_grid()

ax.axis('scaled')
ax.axis([500, 800, 200, 400])

Note that you can pass ax into PlotMapView. Would that work?

rubencalje commented 8 months ago

Hi @langevin-usgs, thank you for the reply.

You are right, I forgot to supply the ax in PlotMapView in my example code.

In my workflow (and possibly of other people as well) I generate a figure with a specified extent to plot a lot of stuff, and also some properties of the model. Then it is useful that the axis limits are not changed by PlotMapView, just like GeoPandas does not change fixed axes limits when a GeoDataFrame is plotted.

This is not a mayor issue, as it is easy to work around, as you show. But it may be a small improvement to the plotting features of flopy.

jlarsen-usgs commented 8 months ago

Hi @rubencalje,

PlotMapView has an exent parameter that can be used to set the plotting extent. The current behaviour, if the extent parameter is not passed in by the user, is to calculate the modelgrid extent (based on coordinates and rotation) and set the axis limits based on that. The reason we explicitly do this is that matplotlib buffers around the data when no axis limits are set.

Here's an example of using the extent parameter:

f, ax = plt.subplots()
ax.axis('scaled')
pmv = flopy.plot.PlotMapView(modelgrid=modelgrid, ax=ax, extent=[500, 800, 200, 400])
pmv.plot_grid()
rubencalje commented 8 months ago

Hi @jlarsen-usgs,

Thanks for your tip on adding the extent parameter to PlotMapView. That is also a good workaround.

My point was that it is possible to check if the user has already set the axes limits, so PlotMapView does not change the axes limits if they are fixed already. I believe you do not think this is worth a change, so let's keep it the way it is now.

Thank you for your time. I will close this issue now.

langevin-usgs commented 8 months ago

Hey @rubencalje, what if we use a special setting for extent that will do what you want. Maybe if the extent is an empty list, for example, we could not change the axes limits in the plot routine? Would that solve your use case?

rubencalje commented 8 months ago

Hey @rubencalje, what if we use a special setting for extent that will do what you want. Maybe if the extent is an empty list, for example, we could not change the axes limits in the plot routine? Would that solve your use case?

That would help.

One last try: a cleaner way would be to check if the axis limits have not been set, with the line ax.get_autoscale_on(), and only set the axes limits to the model extent when this is True.

jlarsen-usgs commented 8 months ago

@rubencalje

I like that suggestion and didn't realize that matplotlib axes objects have a method to check whether they have been autoscaled or previously scaled by the user. I'll put together an update for the plotting routines to incorporate this as a check before we rescale the axes limits.

rubencalje commented 8 months ago

That is great!

If you have little time I can also add a small pull request as well, with some minor changes that should work. I first thought the code could be much simpler, by only adjusting the axes limits in the __init__ method of PLotMapView, if necessary. But in all the plot-methods in PlotMapView it is possible to add a (different) axes as an argument, which makes it slightly more difficult.

I changed

if extent is not None:
    self._extent = extent
else:
    self._extent = None

by

if extent is None:
    if self.ax.get_autoscale_on():
        self._extent = self.mg.extent
    else:
        self._extent = None
else:
    self._extent = extent

and

@property
def extent(self):
    if self._extent is None:
        self._extent = self.mg.extent
    return self._extent

by

@property
def extent(self):
    if self._extent is None:
        return self.ax.axis()
    else:
        return self._extent

Then, in each of the plot methods I replaced:

ax.set_xlim(self.extent[0], self.extent[1])
ax.set_ylim(self.extent[2], self.extent[3])

by

if ax.get_autoscale_on():
    ax.axis(self.extent)
jlarsen-usgs commented 8 months ago

@rubencalje

I drafted up some code last night that I think will solve the issue. Instead of setting axes limits directly in each plotting method, we can route the limit settings through a generalized function like this:

    def _set_axes_limits(self, ax):
        """
        Internal method to set axes limits

        Parameters
        ----------
        ax : matplotlib.pyplot axis
            The plot axis

        Returns
        -------
        ax : matplotlib.pyplot axis object

        """
        if ax.get_autoscale_on():
            ax.set_xlim(self.extent[0], self.extent[1])
            ax.set_ylim(self.extent[2], self.extent[3])
        return ax

I still need to put some tests together to make sure it's working as expected, but I think this issue should be resolved fairly soon.

rubencalje commented 8 months ago

Great! Thank you.