satpy icon indicating copy to clipboard operation
satpy copied to clipboard

Composites combining HRFI and FDHSI channels together

Open guidocioni opened this issue 1 year ago • 28 comments

Hey guys. Today we were trying to create some composites using channels both from HRFI and FDHSI. In this case we were focusing on the ash composite.

The first question I have: where is the definition of the ash composite located? I tried to look into /etc/composites/fci.yaml and /etc/composites/seviri.yaml but couldn't find it. There doesn't seem to be a generic or common composite YAML file, so I'm not sure where they're pulled from (because they are indeed available for both fci and seviri readers). I was interested specifically in knowing which channels are used when composing on FCI, as the channels do not have exactly the same wavelength.

The second question/observation is there may be some error when using channels from the High resolution and Normal resolution together. In this image I'm comparing the ash composite created with channels only from FDHSI (left) and with channels coming from both HRFI and FDHSI. To generate this second image I just provided all files in the scene inputs.

ash

Are these artefacts we're seeing coming from the fact that the data is not yet operational? Or is there some correction to be done when combining data coming from HRFI and FDHSI? I'm not sure which channels are actually used to produce the right image, but there's definitely something from HRFI used.

guidocioni avatar Nov 14 '24 09:11 guidocioni

The definition of ash is here: https://github.com/pytroll/satpy/blob/main/satpy/etc/composites/visir.yaml#L125

I'll post you my "resolution-aware" configs after lunch. There I have set the channel differences to use channels that have equal resolution.

pnuu avatar Nov 14 '24 09:11 pnuu

where is the definition of the ash composite located?

composites/visir.yaml

Are these artefacts we're seeing coming from the fact that the data is not yet operational?

Yes, there are coregistration issues between 10.5 HRFI and other IR channels, so composites that take the difference are affected. Might look better if you use FDHSI only

gerritholl avatar Nov 14 '24 09:11 gerritholl

where is the definition of the ash composite located?

composites/visir.yaml

Are these artefacts we're seeing coming from the fact that the data is not yet operational?

Yes, there are coregistration issues between 10.5 HRFI and other IR channels, so composites that take the difference are affected. Might look better if you use FDHSI only

Ok, I was wondering if this was already known. Well then I guess the only option is to use Normal Resolution data only.

guidocioni avatar Nov 14 '24 10:11 guidocioni

Even with the same resolutions, there are sometimes some weird artifacts (here above Bretagne for example) where it seems that the clouds don't exactly overlap in the different channels. Any idea on how to deal with that?

Reyrem avatar Nov 14 '24 10:11 Reyrem

Even with the same resolutions, there are sometimes some weird artifacts (here above Bretagne for example) where it seems that the clouds don't exactly overlap in the different channels. Any idea on how to deal with that?

That's the co-registration issue Gerrit mentioned.

Here's my ash composite definition, which I've placed in $SATPY_CONFIG_PATH/composites/fci.yaml:

sensor_name: visir/fci

composites:

  ash:
    compositor: !!python/name:satpy.composites.GenericCompositor
    prerequisites:
    - compositor: !!python/name:satpy.composites.DifferenceCompositor
      prerequisites:
      - name: ir_123
        resolution: 2000
      - name: ir_105
        resolution: 2000
    - compositor: !!python/name:satpy.composites.DifferenceCompositor
      prerequisites:
      - name: ir_105
        resolution: 2000
      - name: ir_87
        resolution: 2000
    - name: ir_105
      resolution: 1000
    standard_name: ash

I've done similar changes to other composites that use DifferenceCompositor.

pnuu avatar Nov 14 '24 10:11 pnuu

Even with the same resolutions, there are sometimes some weird artifacts (here above Bretagne for example) where it seems that the clouds don't exactly overlap in the different channels. Any idea on how to deal with that?

That's the co-registration issue Gerrit mentioned.

Here's my ash composite definition, which I've placed in $SATPY_CONFIG_PATH/composites/fci.yaml:

sensor_name: visir/fci

composites:

  ash:
    compositor: !!python/name:satpy.composites.GenericCompositor
    prerequisites:
    - compositor: !!python/name:satpy.composites.DifferenceCompositor
      prerequisites:
      - name: ir_123
        resolution: 2000
      - name: ir_105
        resolution: 2000
    - compositor: !!python/name:satpy.composites.DifferenceCompositor
      prerequisites:
      - name: ir_105
        resolution: 2000
      - name: ir_87
        resolution: 2000
    - name: ir_105
      resolution: 1000
    standard_name: ash

I've done similar changes to other composites that use DifferenceCompositor.

As far as I can see this produces exactly the same image as just using the normal resolution channels without any custom composite. It kind of makes sense, as the resolution for the first 2 bands is forced to be 2km, which is the same as the normal resolution channel. However I would expect to see some differences as the last band is taken at the original 1 km resolution. Probably that's not enough to see a difference in the image.... If the result is the same I would argue it doesn't make sense to get into the trouble of custom defining a composite as just producing the composite with just normal resolution data would give the same product.

guidocioni avatar Nov 14 '24 10:11 guidocioni

There are two effects overlapping here causing the artefacts:

  1. there is indeed an issue in the data that affects the coregistration (overlap) of different IR channels, in this specific case the ir_105-ir_85 pair most predominantly. This is causing most of the colorful "shadow effects" in the images, both when using FDHSI only as well as FDHSI+HRFI. A fix/clear improvement is on the way in upcoming processing facility updates.
  2. there is a further intrinsic problem related to using channels of different resolutions in the same composite, which is particularly evident when doing channel differences (as in the dust/ash RGBs for example). The fact that the HRFI channels contain finer spatial structures that are missing in the FDHSI channels will inevitably lead to artefacts in the images, particularly when doing channel differences, particularly around edgy features. I created a simplified example below trying to explain this in a schematic way - hope it's somewhat clear :D image

ameraner avatar Nov 14 '24 11:11 ameraner

Probably dumb question...would using a different resampling algorithm that adopts more smoothing (linear, bicubic...) alleviate this issue but still keep some spatial features of the high resolution channel?

guidocioni avatar Nov 14 '24 11:11 guidocioni

Ah yes, I forgot to mention that - indeed, a different resampling would likely improve the results a bit, but iirc we tried in the past and it wasn't so much better in the end, still worth investigating though. I guess a very smart resampling algorithm, maybe a machine learning model that uses the HRFI spatial features to enhance the upsampled FDHSI channels, would deliver better results...

ameraner avatar Nov 14 '24 11:11 ameraner

Ah yes, I forgot to mention that - indeed, a different resampling would likely improve the results a bit, but iirc we tried in the past and it wasn't so much better in the end, still worth investigating though. I guess a very smart resampling algorithm, maybe a machine learning model that uses the HRFI spatial features to enhance the upsampled FDHSI channels, would deliver better results...

You heard @Reyrem , now you're next :P

guidocioni avatar Nov 14 '24 11:11 guidocioni

However I would expect to see some differences as the last band is taken at the original 1 km resolution. Probably that's not enough to see a difference in the image....

Yeah, I was also kinda suprised not to see any change between "all 2 km" and "non-difference data at high-res" versions.

If the result is the same I would argue it doesn't make sense to get into the trouble of custom defining a composite as just producing the composite with just normal resolution data would give the same product.

If you want to create any high-res images and include the HRFI files, you do need to have the resolutions defined for at least the channel differences. It might need some more investigation and trials to see whether the single HRFI channel has an effect or not. The effect isn't big, based on what we've seen so far.

Probably dumb question...would using a different resampling algorithm that adopts more smoothing (linear, bicubic...) alleviate this issue but still keep some spatial features of the high resolution channel?

I'm processing the data using bilinear resampling (gradient search) and the results are not any better compared to nearest resampling.

pnuu avatar Nov 14 '24 13:11 pnuu

Ah yes, I forgot to mention that - indeed, a different resampling would likely improve the results a bit, but iirc we tried in the past and it wasn't so much better in the end, still worth investigating though. I guess a very smart resampling algorithm, maybe a machine learning model that uses the HRFI spatial features to enhance the upsampled FDHSI channels, would deliver better results...

This is pretty much what pan sharpening does, isn't it? Technique has been around for a long time and would probably work even better in this case as the spectral band is identical.

simonrp84 avatar Nov 14 '24 13:11 simonrp84

That (pan sharpening) is something I have also been thinking about from time to time and would be a nice addition I think.

BENR0 avatar Nov 19 '24 13:11 BENR0

Here's my ash composite definition, which I've placed in $SATPY_CONFIG_PATH/composites/fci.yaml:

sensor_name: visir/fci

composites:

  ash:
    compositor: !!python/name:satpy.composites.GenericCompositor
    prerequisites:
    - compositor: !!python/name:satpy.composites.DifferenceCompositor
      prerequisites:
      - name: ir_123
        resolution: 2000
      - name: ir_105
        resolution: 2000
    - compositor: !!python/name:satpy.composites.DifferenceCompositor
      prerequisites:
      - name: ir_105
        resolution: 2000
      - name: ir_87
        resolution: 2000
    - name: ir_105
      resolution: 1000
    standard_name: ash

I've done similar changes to other composites that use DifferenceCompositor.

What if HRFI data are missing and you need to produce with FDHSI alone? This will fail, as it will not find ir_105 at 1000 metre resolution:

[DEBUG: 2024-12-02 17:50:51 : satpy.scene] Missing prerequisite for 'DataID(name='dust')': 'DataID(name='ir_105', wavelength=WavelengthRange(min=9.8, central=10.5, max=11.2, unit='µm'), resolution=1000, calibration=<2>, modifiers=())'
[DEBUG: 2024-12-02 17:50:51 : satpy.scene] Missing prerequisite for 'DataID(name='ash')': 'DataID(name='ir_105', wavelength=WavelengthRange(min=9.8, central=10.5, max=11.2, unit='µm'), resolution=1000, calibration=<2>, modifiers=())'

gerritholl avatar Dec 02 '24 16:12 gerritholl

Since I posted that, I actully removed the resolution: 1000 from the green channel. So if HRFI is completely missing, I'd still get an image.

pnuu avatar Dec 02 '24 16:12 pnuu

First test with

  night_microphysical:
    compositor: !!python/name:satpy.composites.GenericCompositor
    prerequisites:
    - compositor: !!python/name:satpy.composites.DifferenceCompositor
      prerequisites:
      - name: ir_123
        resolution: 2000
      - name: ir_105
        resolution: 2000
    - compositor: !!python/name:satpy.composites.DifferenceCompositor
      prerequisites:
      - compositor: !!python/name:satpy.composites.Filler
        prerequisites:
        - name: ir_105
          resolution: 1000
        - name: ir_105
          resolution: 2000
      - compositor: !!python/name:satpy.composites.Filler
        prerequisites:
        - name: ir_38
          resolution: 1000
        - name: ir_38
          resolution: 2000
    - compositor: !!python/name:satpy.composites.Filler
      prerequisites:
      - name: ir_105
        resolution: 1000
      - name: ir_105
        resolution: 2000
    standard_name: night_microphysical

seems promising.

With full set of FDHSI and HRFI files: 20241010_0510_MTG-I1_EPSG_3035_1km_night_microphysical_full-0

Missing HRFI segment number 38: 20241010_0510_MTG-I1_EPSG_3035_1km_night_microphysical_no_hrfi38-0

pnuu avatar Dec 03 '24 06:12 pnuu

And don't mind that the composite doesn't look like what Night Microphysical should, I seemed to use a day-time scene :see_no_evil:

pnuu avatar Dec 03 '24 06:12 pnuu

Looks good! Any significant performance hit?

gerritholl avatar Dec 03 '24 07:12 gerritholl

And what happens if you have no HRFI data at all? Won't it fail to find ir_105 at resolution 1000 m? Or can the Filler deal with this?

gerritholl avatar Dec 03 '24 08:12 gerritholl

No performance hit that I can see. Compared to the normal version that leaves the gap, the processing time is within 0.1 seconds for this one composite and similarly negligible effect on memory usage.

And indeed, if HRFI data are completely missing this will crash.

pnuu avatar Dec 03 '24 08:12 pnuu

Since I posted that, I actully removed the resolution: 1000 from the green channel. So if HRFI is completely missing, I'd still get an image.

Should we include these recipes with satpy? As far as I can see, affected RGBs are:

  • airmass
  • ash
  • convection
  • dust
  • fog
  • night_fog
  • night_microphysics
  • 24h_microphysics

gerritholl avatar Dec 05 '24 10:12 gerritholl

Hello, for what it's worth I found the best result by using the 1km IR channels in the green channel BT difference as well. So, red is 2km/2km, green is 1km/2km, and blue 1km. Using the 1km ir_105 in the red band gives really annoying artifacts but for the green band it is ok.

So, that's a good compromise for having sharpness but avoiding the big "red blobs" when using the 1km ir_105 channel in the red band.

I show below

R : ir_105 at 1 / G: ir_105 at 1km Image

R : ir_105 at 1 / G: ir_105 at 2km Image

R : ir_105 at 2 / G: ir_105 at 1km Image

I plan to use the third option above ("R : ir_105 at 2 / G: ir_105 at 1km") for upcoming work with the ash RGB.

pdebuyl avatar Mar 31 '25 12:03 pdebuyl

Another example, 24h microphysics, close to the western limb:

Image Image

mraspaud avatar Apr 08 '25 08:04 mraspaud

What settings ?

pdebuyl avatar Apr 08 '25 08:04 pdebuyl

What settings ?

Satpy takes the best resolution by default, so the top image is with mixed resolution. Or did I misunderstand the question?

mraspaud avatar Apr 08 '25 11:04 mraspaud

You answered correctly :-) (in the sense, this is the information I was looking for). For the ASH composite, I found that mixing 1km/2km in the R channel created noticeable artifacts while mixing 1km/2km in the G channel was ok. It will depend composite by composite of course, but I was happy to have "high resolution ash while reducing artifacts". My tests also seemed to indicate that gradient_search produced better-looking images than nearest.

pdebuyl avatar Apr 08 '25 12:04 pdebuyl

You answered correctly :-) (in the sense, this is the information I was looking for). For the ASH composite, I found that mixing 1km/2km in the R channel created noticeable artifacts while mixing 1km/2km in the G channel was ok. It will depend composite by composite of course, but I was happy to have "high resolution ash while reducing artifacts". My tests also seemed to indicate that gradient_search produced better-looking images than nearest.

Hey @pdebuyl would you mind posting the full YAML you ended up using to create the ash composite at the highest possible resolution while reducing the artifacts? I read all the posts before but I'm confused now :)

guidocioni avatar Apr 09 '25 12:04 guidocioni

I use the following $SATPY_CONFIG_PATH/composites/fci.yaml

sensor_name: visir/fci

composites:

  ash_22_12_1:
    compositor: !!python/name:satpy.composites.GenericCompositor
    prerequisites:
    - compositor: !!python/name:satpy.composites.DifferenceCompositor
      prerequisites:
      - name: ir_123
        resolution: 2000
      - name: ir_105
        resolution: 2000
    - compositor: !!python/name:satpy.composites.DifferenceCompositor
      prerequisites:
      - name: ir_105
        resolution: 1000
      - name: ir_87
        resolution: 2000
    - name: ir_105
      resolution: 1000
    standard_name: ash_22_12_1

And this for the $SATPY_CONFIG_PATH/enhancements/fci.yaml

enhancements:
  ash_22_12_1:
    standard_name: ash_22_12_1
    operations:
    - name: stretch
      method: !!python/name:satpy.enhancements.stretch
      kwargs:
        stretch: crude
        min_stretch: [-4, -4, 243]
        max_stretch: [2, 5, 303]

Note: I probably don't need to duplicate this in the enhancement file but I found this solution to make sure to have the proper stretching.

pdebuyl avatar Apr 09 '25 13:04 pdebuyl