GenericMappingTools / pygmt

A Python interface for the Generic Mapping Tools.
https://www.pygmt.org
BSD 3-Clause "New" or "Revised" License
758 stars 220 forks source link

Does a "clear(figure)" instance exists ? #2111

Open Sarenrenegade opened 2 years ago

Sarenrenegade commented 2 years ago

Description of the desired feature

Hello, I'm currently working on a program designed to make an anilation of seismic wave propagation using pygmt : Basically with seismograms on Z component from many stations of an earthquake I have the displacement at every time step. Since there is no "pygmt.animate ()" I decided to plot multiple maps for each steps and then make the animation with opencv. It works great when I have few figures to plot, but for the animation to be smooth (over 1000 images) my computer runs out of memory. My current workaround is to launch the program again from the last step. Is there a way to use plain gmt in pygmt to end the session/close a figure to free memory ? The idea being to have all the process from downloading seismograms to making the animation in python I don't want to use plain gmt in cmd. I also don't want to use cartopy: I find pygmt easier to work with. Thanks

Are you willing to help implement and maintain this feature? Maybe

I'm not really sure I have the skills, I rarely use python the "right" way (with objects not sequentially) and I have never worked on a wrapper , but I'm willing to help

weiji14 commented 2 years ago

Hi @Sarenrenegade, thanks for trying out PyGMT! What you described on creating an animation seems to be the job for https://docs.generic-mapping-tools.org/6.4/movie.html (see also GMT animation gallery at https://docs.generic-mapping-tools.org/6.4/animations.html and feature request at https://github.com/GenericMappingTools/pygmt/issues/1347). Unfortunately, this movie function isn't wrapped in PyGMT yet, but it's something that I'm sure a lot of people would love to see in PyGMT!

In the meantime, there is a way to call unwrapped GMT C modules in PyGMT via https://www.pygmt.org/v0.7.0/api/generated/pygmt.clib.Session.call_module.html. This is slightly better than calling GMT with Python's os.system or subprocess, but not as nice as a fully wrapped PyGMT module. As an example, you could use:

import pygmt

with pygmt.clib.Session() as lib:
    # movie(canvassize="480p", nframes=49, displayrate=3, format="gif", prefix="frame", erase=True)
    lib.call_module("movie", "animation.sh -C480p -T49 -D3 -A+l -Nframe -Z")

Now I've not tested if this works, but it's based on some code I have at https://github.com/weiji14/30DayMapChallenge2021/blob/8dc7fa4597e8b98a908071243a92239c1aec54b8/day20_movement.py#L99-L125 (which you're more than welcome to modify). Only downside is that you'll need to have an animate.sh shell script.

Is there a way to use plain gmt in pygmt to end the session/close a figure to free memory ?

In fact, we do have have a begin and end function, though it's 1) not documented and 2) a little dangerous to use unless you know what you're doing. You can see it being used here:

https://github.com/GenericMappingTools/pygmt/blob/e5f5da9ba4713912ee96bd5eed80aa697d40de46/pygmt/tests/test_session_management.py#L7-L23

Let us know though if you need help with this, it'll be easier if you post a short example code so we know what we're working with.

Oh, and you're also more than welcome to help implement the movie function. We've got some draft instructions on how to wrap a new module at #1687, specifically at https://github.com/GenericMappingTools/pygmt/blob/wrap-module-instructions/doc/contributing.md#initial-feature-implementation :wink:

Sarenrenegade commented 2 years ago

Thanks for your answer. I will try your solution for the animation and keep you aprised. I will also take a look at the instructions to wrap a module (no promises!).

The (simplified) function I use to generate the images is below : Do you think using gmt end in my loop is wise? For now if the program crashes the user has the opportunity to start from the last generated image, but it's not exactly user friendly...

I included a .txt file that you can use as an input, if you want to try.

def plot_map(trace, region_of_interest):

    """ This function will generate one frame for each step point of the seismogram for n stations

    INPUT : array (n_lines, m_columns), column 1 : longitude, column 2: latitude, column 3 to m: seismogram, type : numpy array
                 region of interest : [min longitude, max long, min lat, max lat], type: list

    OUTPUT : n images classified sequentially 

    !!! warning !!! : for the video to work the images must be named sequentially (eg : 0001 for the  first image if you have over 1000 image for the video """ 

    lon=trace[:,0]
    lat=trace[:,1]
    t=trace[:,2:trace.shape[1]-2]

    fig = pygmt.Figure()
    fig.basemap(region=region_of_interest, projection="M15c", frame=True)
    fig.coast(land="grey")
    pygmt.makecpt(cmap="polar", series=[-1.2*1e4, 1.2*1e4]) # Amatrice filtered twice : 1.2*1e4

      for n in range (0, t.shape[1]):
          fig.plot(x=lon, y=lat,color=t[:,n],cmap=True, style="c0.1c")
          fig.colorbar(position="JMR+o0.5c/0c+w5c/0.3c+n", scale=1e-3, frame=['af+ldéplacement_(z)','y+lmm/s2'])

          if n<10:
              fig.savefig('filtered_00000'+str(n)+'.png')
          elif 10<=n<=99:
              fig.savefig('filtered_0000'+str(n)+'.png')
          elif 100<=n<=999:
              fig.savefig('filtered_000'+str(n)+'.png')
          elif 1000<=n<=9999:
              fig.savefig('filtered_00'+str(n)+'.png')
          elif 10000<=n<=99999:
              fig.savefig('filtered_0'+str(n)+'.png')
          elif 100000<=n<=999999:
              fig.savefig('filtered_'+str(n)+'.png')

import numpy as np
trace= np.loadtxt('Amatrice_2016_V2__filtered_normalized.txt'
region=[-11, 30, 35, 56]
plot_map(trace, region)

Amatrice_2016_V2__filtered_normalized.txt

weiji14 commented 2 years ago

The (simplified) function I use to generate the images is below : Do you think using gmt end in my loop is wise?

If the for-loop is run one-by-one, then calling gmt end should be ok (but you'll need to be careful to call gmt begin at the start). To be safe, you might want to have something like this (as mentioned in https://github.com/GenericMappingTools/pygmt/issues/2111#issuecomment-1238169752):

import pygmt
from pygmt.session_management import begin, end 

def plot_map():
    end()
    begin()

    fig = pygmt.FIgure()
    # Plotting code here

    end()
    begin()

Now if you're trying to run the function in parallel, see https://github.com/GenericMappingTools/pygmt/issues/217#issuecomment-754774875.

To be honest, I would actually recommend you to use gmt movie if possible, because that should be the most stable in terms of memory management (less crashes) and the -x option would allow you to create each frame in parallel in a safe and fast way. However, it might not have the 'restart from last image' feature like your code above :slightly_smiling_face:

P.S. You can use str(n).zfill(6) to pad your filename with leading zeros. See https://stackoverflow.com/questions/339007/how-do-i-pad-a-string-with-zeroes.

Sarenrenegade commented 2 years ago

If the for-loop is run one-by-one

Yes it is, since the program is mainly designed for education I wasn't exactly looking for efficiency. But it's good to know that there is a faster solution. Thanks for sharing! And I'm pretty sure I can mke it work with the start from last image feature ;-)

Also thanks for the tip for the leading zeros!!

Anyway, I will definetly read the doc on how to add a module in pygmt, because it would be really great to have an animate module.