Increasing speed/parallel processing of Satpy save_datasets

[Fabio Mingione]

Dear Pytroll Community,

I would like to make my Satpy MSG 0° processing chain faster.

In particular, we're looking for an internal Satpy command to make the save_datasets method run faster.

Looking at the Satpy documentation pdf, pg. 361, we tried the num_threads=4 option (default should be 1 as far as I understood), but nothing changed in terms of execution time. I saw num_threads was referred to GeoTIFFs in particular, should it work with MSG as well? Is there any similar command to accomplish our goal? 

Furthermore, since HRIT files arrive at our machines every 15mins and we launch several python scripts (saving different scenes) at the same time, I noticed a delay is accumulated at each time slot, causing the chain to miss a slot once in a while. We partially improved the speed by running the scripts in parallel rather than in series, so the internal Satpy command --in case it exists-- would only be an addition to this. 

Thanks a lot!

Regards,

Fabio

[lobsiger...@gmail.com]

Fabio,

here is my 2 cents for speeding up MSG4 processing:

- Depending on the different composites you are using you should only  decompress the channel files you need and only decompress them once.

- There are two possibilities of caching GEO sat processing chains:

  a) There is a area depending "projection cache" when resampling:
  new_scene = global_scene.resample(area,  radius_of_influence= 20000,  cache_dir = some_place)
  b) There is an overlay cache for coastlines/grid etc. (only useful if pycoast overlays are used at all).
  new_scene.save_dataset( .....'overlays': {...., 'cache': {'file': ..., 'regenerate': ...}})

Regards,
Ernst

[David Hoese]

Fabio,

It wasn't clear to me from your message what format you are saving to?
The reading and the writing in satpy are two separate things so writing
to a geotiff and reading from a MSG file are relatively separate tasks
as far as what options are available to control them. So some questions
and other thoughts to consider:

1. What format are you writing to?
2. If geotiff, are you tiling the geotiffs?
3. How many CPU cores and how much memory does your system have?
4. How many dask workers are you using?
5. Are you reading from or writing to network storage?
6. Are you doing any other optimizations or environment variable setting
to try to speed things up (ex. PYTROLL_CHUNK_SIZE)?
6. I agree with Ernst on the resample caching.

As for overall geotiff performance, we are using rasterio to write the
geotiffs and it is relatively fast if you specify the "right" options
and have well chunked dask array operations. We will in the future
switch to using rioxarray for writing geotiffs which has a little more
flexibility with how geotiffs are written on a distributed/cluster
system, but for single machine processing I don't think this should make
a difference.

Ernst,

Have you see a lot of improvement when you use the overlay caching? In
my own very basic tests I don't remember it saving me more than half a
second or so compared to running it without caching.

Dave

> --
> You received this message because you are subscribed to the Google
> Groups "pytroll" group.
> To unsubscribe from this group and stop receiving emails from it, send
> an email to pytroll+u...@googlegroups.com
> <mailto:pytroll+u...@googlegroups.com>.
> To view this discussion on the web, visit
> https://groups.google.com/d/msgid/pytroll/6399f9cb-996d-4341-8d92-9e9bcb4d73ddn%40googlegroups.com
> <https://groups.google.com/d/msgid/pytroll/6399f9cb-996d-4341-8d92-9e9bcb4d73ddn%40googlegroups.com?utm_medium=email&utm_source=footer>.

[Fabio Mingione]

Ernst, Dave, thank you very much for your replies! 

I will conduct some tests in the next days based on these.

In the meanwhile, @Dave:

Sorry for the lack of clarity:

1. I'm writing to JPG

2. I'm not using geotiffs. I'm reading HRIT and writing JPG.

3. 16 CPU cores, 32GB RAM.

4. I never modified DASK_NUM_WORKERS. I read in the documentation that it's suggested to set it between 4 and 8, right?

5. I'm reading and writing locally. All necessary HRITs arrive already decompressed at the Satpy input location and Satpy starts processing only once that is completed. 

6. No. I never modified PYTROLL_CHUNK_SIZE, so it defaults to 4096. I will try to set it to 2048 and increase DASK_NUM_WORKERS.

7. I will try the resample caching too.

I will get back you if I find anycombination of these to decrease the total computation time.

Thanks,

Fabio

[David Hoese]

Fabio,

1. A single JPEG file can't be written in parallel (same for all simple
image formats). There isn't a lot you'll be able to do to speed up the
writing of these.
4. 4 and 8 dask workers have worked well for a lot of us.
6. Dask will default to the number of cores you have on your system. In
some cases this can create too many workers and make it difficult for
your system. You may want to look at some of the suggestion in the satpy
FAQ, specifically
https://satpy.readthedocs.io/en/latest/faq.html#why-multiple-cpus-are-used-even-with-one-worker
7. Note that the resampling cache does not do anything if you are doing
"native" resampling.

Dave

> <https://groups.google.com/d/msgid/pytroll/6399f9cb-996d-4341-8d92-9e9bcb4d73ddn%40googlegroups.com?utm_medium=email&utm_source=footer
> <https://groups.google.com/d/msgid/pytroll/6399f9cb-996d-4341-8d92-9e9bcb4d73ddn%40googlegroups.com?utm_medium=email&utm_source=footer>>.

>
>
> --
> You received this message because you are subscribed to the Google
> Groups "pytroll" group.
> To unsubscribe from this group and stop receiving emails from it, send
> an email to pytroll+u...@googlegroups.com
> <mailto:pytroll+u...@googlegroups.com>.
> To view this discussion on the web, visit

> https://groups.google.com/d/msgid/pytroll/e39c0025-3d10-440f-a5cb-88953931636en%40googlegroups.com
> <https://groups.google.com/d/msgid/pytroll/e39c0025-3d10-440f-a5cb-88953931636en%40googlegroups.com?utm_medium=email&utm_source=footer>.

[Fabio Mingione]

Thank you Dave! 

Now I understand a few things better

[lobsiger...@gmail.com]

Fabio,

reorganizing my scripts I did some speed testing. As you said
you make different composites, one point that -- as far as I
can see -- has not yet been discussed above. TBH this was rather
new to me as well: If you produce different composites that
you resample to the same area you can load a full composite list
and delay composite generation until you save your images. This
way the different channel datasets are only resampled once which
is time saving if you have a list of composits that share a number
of channels. Below I resampled MSG4 data to my area ='westminster'
(stere, 2048 x 2048 pixels) using the following composite list:

composites = ['natural_color', 'dust', 'fog', 'green_snow', 'ir_sandwich']

I made two scripts that work under Devuan Beowulf GNU/Linux as follows:

MSG4-S.py

for composite in composites:
    files = find all files necessary for the single composite at hand
    global_scene = Scene(filenames = files, reader = 'seviri_l1b_hrit')
    global_scene.load([composite], upper_right_corner = 'NE')
    new_scene = global_scene.resample(area, radius_of_influence = 20000)
    my_overlays = produce overlay dict for coastlines, rivers, borders, points, ...
    new_scene.save_dataset(composite, 'MSG4.png', overlay = {'coast_dir': dbroot, 'overlays' : my_overlays})
    Additional ImageMagick step to produce legend at left and save final result as *.jpg


MSG4-X.py

files = find all files necessary for the whole list of composites
global_scene = Scene(filenames = files, reader = 'seviri_l1b_hrit')
global_scene.load(composites, upper_right_corner = 'NE', generate = False)
new_scene = global_scene.resample(area, radius_of_influence = 20000)
my_overlays = produce overlay dict for coastlines, rivers, borders, points, ...
for composite in composites:
    new_scene.save_dataset(composite, 'MSG4.png', overlay = {'coast_dir': dbroot, 'overlays' : my_overlays})
    Additional ImageMagick step to produce legend at left and save final result as *.jpg


I get the following results on a 10 year old I7 with 16GB of RAM that is
used as an EUMETCast receiver and also makes images and movies croned:

MSG4-S.py

1) All files used by every composite in the loop have to be decompressed
   This way some/many(?) files have to be xRITDecompressed several times
   SatPy caching: Proj-cache=False (GEO sats only), Overlay-cache=False
   Time (real) for generation of all 5 composites: 1Min 37sec

2) Decompressed files are cached if already used by some composite before
   This way files needed in the loop have to be xRITDecompressed only once
   SatPy caching: Proj-cache=False (GEO sats only), Overlay-cache=False
   Time (real) for generation of all 5 composites: 1Min 26sec

3) All xrit files needed in the loop have been decompressed beforehand
   SatPy caching: Proj-cache=False (GEO sats only), Overlay-cache=False
   Time (real) for generation of all 5 composites: 1Min 11sec


MSG4-X.py

4) All xrit files needed for composites have been decompressed beforehand
   SatPy caching: Proj-cache=False (GEO sats only), Overlay-cache=False
   Time (real) for generation of all 5 composites: 56sec

5) All xrit files needed for composites have been decompressed beforehand
   SatPy caching: Proj-cache=True (GEO sats only), Overlay-cache=False
   Time (real) for generation of all 5 composites: 40sec

6) All xrit files needed for composites have been decompressed beforehand
   SatPy caching: Proj-cache=True (GEO sats only), Overlay-cache=True
   Time (real) for generation of all 5 composites: 35sec


Results are attached below. Original image size is 2401 x 2048 pixels.
Maybe MSG4-X.py is what you already do. For me it's an big step forward.


Best regards,
Ernst

[Panu Lahtinen]

On 20.6.2022 15.01, lobsiger...@gmail.com wrote:
>
> MSG4-X.py
>
> files = find all files necessary for the whole list of composites
> global_scene = Scene(filenames = files, reader = 'seviri_l1b_hrit')
> global_scene.load(composites, upper_right_corner = 'NE', generate = False)
> new_scene = global_scene.resample(area, radius_of_influence = 20000)
> my_overlays = produce overlay dict for coastlines, rivers, borders,
> points, ...
> for composite in composites:
>     new_scene.save_dataset(composite, 'MSG4.png', overlay =
> {'coast_dir': dbroot, 'overlays' : my_overlays})
>     Additional ImageMagick step to produce legend at left and save
> final result as *.jpg

Saving is faster without looping the composites. That way the re-used
channels are processed only once.


Panu, a core-dev

--
Finnish Meteorological Institute
Erik Palménin aukio 1, P.O.Box 503, FIN-00101 Helsinki, Finland
Room: 4B15a, Phone: +358 50 380 3261

[lobsiger...@gmail.com]

Panu,

I already see a remarkable speedup in my MSG4-X.py version. You tell me this works even faster?

How do I save the different images without a loop? Does your post mean I can even do something like:

new_scene.save_dataset(composites, list_of_filenames, overlay =  {'coast_dir': dbroot,  'overlays' : my_overlays} )

??????

If so, is there also a way to add different overlays without a loop?

Regards,

Ernst

[Panu Lahtinen]

On 20.6.2022 15.37, lobsiger...@gmail.com wrote:
> Panu,
>
> I already see a remarkable speedup in my MSG4-X.py version. You tell me
> this works even faster?
> How do I save the different images without a loop? Does your post mean I
> can even do something like:
>
> new_scene.save_dataset(composites, list_of_filenames, overlay =
> {'coast_dir': dbroot,  'overlays' : my_overlays} )

You can either rely on the built-in filename pattern:

new_scene.save_datasets(base_dir='/tmp')

or define it yourself:

new_scene.save_datasets(base_dir='/tmp',
filename="{start_time:%Y%m%d_%H%M}_{platform_name}.tif")


P

[lobsiger...@gmail.com]

Panu,

thanks for the quick reply. Instead of .save_dataset() I should use .save_datasets()

One little 's' makes the difference. There remains the question: How can I handle

composite dependant overlays with that. And I will certainly need my loop for

image annotation using ImageMagick. I'll give it a try and report back.

Ernst

[Panu Lahtinen]

On 20.6.2022 16.07, lobsiger...@gmail.com wrote:
> Panu,
>
> thanks for the quick reply. Instead of .save_dataset() I should use
> .save_datasets()
> One little 's' makes the difference. There remains the question: How can
> I handle
> composite dependant overlays with that. And I will certainly need my
> loop for
> image annotation using ImageMagick. I'll give it a try and report back.

Oh, right.

Another thing to try, that might also work with overlays:


from satpy.writers import compute_writer_results

pattern = "{start_time:%Y%m%d_%H%M}_{platform_name}.tif"

save_objects = []
for composite in composites:
overlay_config = ...
save_objects.append(
new_scene.save_dataset(composite,
base_dir='/tmp',
filename=pattern,
overlay=overlay_config,
compute=False)
compute_writer_results(save_objects)

[lobsiger...@gmail.com]

Panu,

I tried the multi save method like:

new_scene.save_datasets(basedir='.', filename='{standard_name}.png',

    overlay = {'coast_dir': dbroot, 'overlays' : my_overlays})

and added a last loop over composites to make the annotation with ImageMagick.

The (real) time I use now for all 5 above attached composites is:       21sec (!!)

I see no possibility to make my_overlays dependant on {standard_name} of

the dataset that is currently saved. So I'll maybe have a look at your proposal

above. You are pushing me to my mental limits ;-). Many thanks for your help.

Ernst

[lobsiger...@gmail.com]

Panu,
Developers,

@Panu I have your save_objects method working with 21s for the
five 'westminster' composites as before. Many thanks for this!

What I see is that composites with BlackMarble backgrounds
do not work with generate=False. This is probably common
to all composites that use background image files. Is there
any chance to workaround this issue in the following ways?

a) Prepare the background file by resampling it to the layout
   of the input data. For MSG4 I would e.g. prepare a BlackDisk
   file 3712x3712 pixels (seviri 'native') by resampling the
   BlackMarble file from 'Plate Carée' to 'GEOS 0deg'. This
   file should stay on disk and work for all resampling areas.

b) Give the simple image PIL writer a first image to put the
   resampled satellite area image on. For my case described
   above with area 'westminster' I would prepare once for all

   and write to disk a BlackWestminster.png image beforehand.

c) Have I missed something and BlackMarbles should work?

Any other ideas for speedup composites with backgrounds?

Best regards,
Ernst

[lobsiger...@gmail.com]

Dave and All,

<cite>
Ernst,

Have you seen a lot of improvement when you use the overlay caching? In

my own very basic tests I don't remember it saving me more than half a
second or so compared to running it without caching.

Dave

</cite>

To answer above question and increase my understanding I have made a
couple of tests. I have put my findings in the text file below. Under the
following wetransfer link (valid one week) you find the text with all
resulting pictures. The bottom line is that overlay caching is also
effective especially if you use coast lines in high resolution areas.
For quality comparison image names bear the resampler used as suffix.

https://we.tl/t-IrVSiFjt2Z

(32 MB *.zip file)

Best regards,
Ernst