Confusing (wrong?) band names for subset with snappy

[Nobody]

Hello.
I don’t find anything like S1_radiance_an in the band names that snappy reads for the “subset” operator.
What I’m doing is reading the .xml of the following S3B_SL_1_RBT____20190219T092217_20190219T092517_20190219T113119_0179_022_150_2160_MAR_O_NR_003.SEN3, resampling everyhting to 500m, reprojecting. Now I want to create two subsets with the bands I need and export to Geotiff.
My problem arises when I try to create the subset with the radiances (S1-S6)_radiance_an. The error I get is:
RuntimeError: org.esa.snap.core.gpf.OperatorException: Operator ‘SubsetOp’: Value for ‘Source Bands’ is invalid: ‘S1_radiance_an’
When I print the list of band names snappy sees (from my reprojected product) it gives me the following, which is quite different from the bands SNAP gives me:
band 0 = x_in
band 1 = y_in
band 2 = x_io
band 3 = y_io
band 4 = F1_BT_in
band 5 = F1_exception_in
band 6 = F1_BT_io
band 7 = F1_exception_io
band 8 = F2_BT_in
band 9 = F2_exception_in
band 10 = F2_BT_io
band 11 = F2_exception_io
band 12 = bayes_in
band 13 = cloud_in
band 14 = confidence_in
band 15 = pointing_in
band 16 = probability_cloud_dual_in
band 17 = probability_cloud_single_in
band 18 = bayes_io
band 19 = cloud_io
band 20 = confidence_io
band 21 = pointing_io
band 22 = probability_cloud_dual_io
band 23 = probability_cloud_single_io
band 24 = elevation_in
band 25 = latitude_in
band 26 = longitude_in
band 27 = elevation_io
band 28 = latitude_io
band 29 = longitude_io
band 30 = detector_in
band 31 = pixel_in
band 32 = scan_in
band 33 = detector_io
band 34 = pixel_io
band 35 = scan_io
band 36 = S7_BT_in
band 37 = S7_exception_in
band 38 = S7_BT_io
band 39 = S7_exception_io
band 40 = S8_BT_in
band 41 = S8_exception_in
band 42 = S8_BT_io
band 43 = S8_exception_io
band 44 = S9_BT_in
band 45 = S9_exception_in
band 46 = S9_BT_io
band 47 = S9_exception_io
band 48 = x_tx
band 49 = y_tx
band 50 = latitude_tx
band 51 = longitude_tx
band 52 = sat_azimuth_tn
band 53 = sat_path_tn
band 54 = sat_zenith_tn
band 55 = solar_azimuth_tn
band 56 = solar_path_tn
band 57 = solar_zenith_tn
band 58 = sat_azimuth_to
band 59 = sat_path_to
band 60 = sat_zenith_to
band 61 = solar_azimuth_to
band 62 = solar_path_to
band 63 = solar_zenith_to
band 64 = cloud_fraction_tx
band 65 = dew_point_tx
band 66 = east_west_stress_tx_time_1_tx
band 67 = east_west_stress_tx_time_2_tx
band 68 = east_west_stress_tx_time_3_tx
band 69 = east_west_stress_tx_time_4_tx
band 70 = east_west_stress_tx_time_5_tx
band 71 = latent_heat_tx_time_1_tx
band 72 = latent_heat_tx_time_2_tx
band 73 = latent_heat_tx_time_3_tx
band 74 = latent_heat_tx_time_4_tx
band 75 = latent_heat_tx_time_5_tx
band 76 = north_south_stress_tx_time_1_tx
band 77 = north_south_stress_tx_time_2_tx
band 78 = north_south_stress_tx_time_3_tx
band 79 = north_south_stress_tx_time_4_tx
band 80 = north_south_stress_tx_time_5_tx
band 81 = sea_ice_fraction_tx
band 82 = sea_surface_temperature_tx
band 83 = sensible_heat_tx_time_1_tx
band 84 = sensible_heat_tx_time_2_tx
band 85 = sensible_heat_tx_time_3_tx
band 86 = sensible_heat_tx_time_4_tx
band 87 = sensible_heat_tx_time_5_tx
band 88 = skin_temperature_tx
band 89 = snow_albedo_tx
band 90 = snow_depth_tx
band 91 = soil_wetness_tx
band 92 = solar_radiation_tx_time_1_tx
band 93 = solar_radiation_tx_time_2_tx
band 94 = solar_radiation_tx_time_3_tx
band 95 = solar_radiation_tx_time_4_tx
band 96 = solar_radiation_tx_time_5_tx
band 97 = specific_humidity_tx_pressure_level_1_tx
band 98 = specific_humidity_tx_pressure_level_2_tx
band 99 = specific_humidity_tx_pressure_level_3_tx
band 100 = specific_humidity_tx_pressure_level_4_tx
band 101 = specific_humidity_tx_pressure_level_5_tx
band 102 = specific_humidity_tx_pressure_level_6_tx
band 103 = specific_humidity_tx_pressure_level_7_tx
band 104 = specific_humidity_tx_pressure_level_8_tx
band 105 = specific_humidity_tx_pressure_level_9_tx
band 106 = specific_humidity_tx_pressure_level_10_tx
band 107 = specific_humidity_tx_pressure_level_11_tx
band 108 = specific_humidity_tx_pressure_level_12_tx
band 109 = specific_humidity_tx_pressure_level_13_tx
band 110 = specific_humidity_tx_pressure_level_14_tx
band 111 = specific_humidity_tx_pressure_level_15_tx
band 112 = specific_humidity_tx_pressure_level_16_tx
band 113 = specific_humidity_tx_pressure_level_17_tx
band 114 = specific_humidity_tx_pressure_level_18_tx
band 115 = specific_humidity_tx_pressure_level_19_tx
band 116 = specific_humidity_tx_pressure_level_20_tx
band 117 = specific_humidity_tx_pressure_level_21_tx
band 118 = specific_humidity_tx_pressure_level_22_tx
band 119 = specific_humidity_tx_pressure_level_23_tx
band 120 = specific_humidity_tx_pressure_level_24_tx
band 121 = specific_humidity_tx_pressure_level_25_tx
band 122 = surface_pressure_tx
band 123 = temperature_profile_tx_pressure_level_1_tx
band 124 = temperature_profile_tx_pressure_level_2_tx
band 125 = temperature_profile_tx_pressure_level_3_tx
band 126 = temperature_profile_tx_pressure_level_4_tx
band 127 = temperature_profile_tx_pressure_level_5_tx
band 128 = temperature_profile_tx_pressure_level_6_tx
band 129 = temperature_profile_tx_pressure_level_7_tx
band 130 = temperature_profile_tx_pressure_level_8_tx
band 131 = temperature_profile_tx_pressure_level_9_tx
band 132 = temperature_profile_tx_pressure_level_10_tx
band 133 = temperature_profile_tx_pressure_level_11_tx
band 134 = temperature_profile_tx_pressure_level_12_tx
band 135 = temperature_profile_tx_pressure_level_13_tx
band 136 = temperature_profile_tx_pressure_level_14_tx
band 137 = temperature_profile_tx_pressure_level_15_tx
band 138 = temperature_profile_tx_pressure_level_16_tx
band 139 = temperature_profile_tx_pressure_level_17_tx
band 140 = temperature_profile_tx_pressure_level_18_tx
band 141 = temperature_profile_tx_pressure_level_19_tx
band 142 = temperature_profile_tx_pressure_level_20_tx
band 143 = temperature_profile_tx_pressure_level_21_tx
band 144 = temperature_profile_tx_pressure_level_22_tx
band 145 = temperature_profile_tx_pressure_level_23_tx
band 146 = temperature_profile_tx_pressure_level_24_tx
band 147 = temperature_profile_tx_pressure_level_25_tx
band 148 = temperature_tx
band 149 = thermal_radiation_tx_time_1_tx
band 150 = thermal_radiation_tx_time_2_tx
band 151 = thermal_radiation_tx_time_3_tx
band 152 = thermal_radiation_tx_time_4_tx
band 153 = thermal_radiation_tx_time_5_tx
band 154 = total_column_ozone_tx
band 155 = total_column_water_vapour_tx
band 156 = u_wind_tx_time_1_tx
band 157 = u_wind_tx_time_2_tx
band 158 = u_wind_tx_time_3_tx
band 159 = u_wind_tx_time_4_tx
band 160 = u_wind_tx_time_5_tx
band 161 = v_wind_tx_time_1_tx
band 162 = v_wind_tx_time_2_tx
band 163 = v_wind_tx_time_3_tx
band 164 = v_wind_tx_time_4_tx
band 165 = v_wind_tx_time_5_tx
I’m not so confident as to which I should choose to export my radiances so, could you please help me? Which ones of the above correspond to the bands S1_radiance_an to S6_radiance_an?
Is it normal to get such names and not the ones I see on SNAP, following the same procedure?

Thank you.

[Nobody]

Ok. I’m not quite sure what the problem was, since I don’t remember having changed anyting regarding the resampling function in my script, however now it gets 285 bands and it does find an S1_radiance_an (and all the rest).
I have the impression, but no proof yet, that sometimes the same functions give me error and some time later they don’t, but I haven’t figured out yet if it’s me getting crazy or it really happens
Anyway, thank you for your time.

[Nobody]

Ok. Now I’m sure. Something is wrong somewhere. The operators I’m ineterested in are Resample, Reproject, Subset and Rad2Ref. So, there are a few arguments I’d like to discuss…

1. I have the following script:

#!/usr/bin/python

coding: utf8

import snappy
from snappy import GPF

‘’’
def get_snap_info(operator):
“”"
Returns information about SNAP operators and their parameters
“”"
op_spi = GPF.getDefaultInstance().getOperatorSpiRegistry().getOperatorSpi(operator)
print(‘Op name:’, op_spi.getOperatorDescriptor().getName())
print(‘Op alias:’, op_spi.getOperatorDescriptor().getAlias())
param_Desc = op_spi.getOperatorDescriptor().getParameterDescriptors()
for param in param_Desc:
print(param.getName(), “or”, param.getAlias())

print “BandMaths”
get_snap_info(“BandMaths”)
‘’’
inp = snappy.ProductIO.readProduct("/MyPath/S3B_SL_1_RBT____20190219T092217_20190219T092517_20190219T113119_0179_022_150_2160_MAR_O_NR_003.SEN3/xfdumanifest.xml")

2. Set parameters

Rad2Refl:

HashMap = snappy.jpy.get_type(‘java.util.HashMap’)
parRad = HashMap()
parRad.put(“sensor”, “SLSTR_500m”)
parRad.put(“conversionMode”, “RAD_TO_REFL”)
parRad.put(“copyTiePointGrids”, “false”)
parRad.put(“copyFlagBandsAndMasks”, “false”)
parRad.put(“copyNonSpectralBands”, “false”)

Resample:

parResam = HashMap()
parResam.put(“targetWidth”, 3000)
parResam.put(“targetHeight”, 2400)
parResam.put(“upsampling”, “Nearest”)
parResam.put(“downsampling”, “First”)
parResam.put(“flagDownsampling”, “First”)
parResam.put(“resampleOnPyramidLevels”, “true”)

Reproject:

parProj = HashMap()
parProj.put(“crs”, ‘PROJCS[“WGS 84 / UTM zone 32N”, GEOGCS[“WGS 84”, DATUM[“World Geodetic System 1984”, SPHEROID[“WGS 84”, 6378137.0, 298.257223563, AUTHORITY[“EPSG”,“7030”]], AUTHORITY[“EPSG”,“6326”]], PRIMEM[“Greenwich”, 0.0, AUTHORITY[“EPSG”,“8901”]], UNIT[“degree”, 0.017453292519943295], AXIS[“Geodetic longitude”, EAST], AXIS[“Geodetic latitude”, NORTH], AUTHORITY[“EPSG”,“4326”]], PROJECTION[“Transverse_Mercator”, AUTHORITY[“EPSG”,“9807”]], PARAMETER[“central_meridian”, 9.0], PARAMETER[“latitude_of_origin”, 0.0], PARAMETER[“scale_factor”, 0.9996], PARAMETER[“false_easting”, 500000.0], PARAMETER[“false_northing”, 0.0], UNIT[“m”, 1.0], AXIS[“Easting”, EAST], AXIS[“Northing”, NORTH], AUTHORITY[“EPSG”,“32632”]]’)
parProj.put(“resampling”, “Nearest”)
parProj.put(“orthorectify”, “false”)
parProj.put(“noDataValue”, -32768.0)
parProj.put(“includeTiePointGrids”, “true”)
parProj.put(“addDeltaBands”, “false”)

Subset:

parSub1 = HashMap()
parSub1.put(“bandNames”, “S1_reflectance_an,S2_reflectance_an,S3_reflectance_an,S5_reflectance_an,S6_reflectance_an”)
parSub1.put(“copyMetadata”, “true”)

outResam = GPF.createProduct(“Resample”, parResam, inp)
outRad = GPF.createProduct(“Rad2Refl”, parRad, outResam)
outSub1 = GPF.createProduct(“Subset”, parSub1, outRad)
outProj = GPF.createProduct(“Reproject”, parProj, outSub1)
outSub2 = GPF.createProduct(“Subset”, parSub1, outProj)

print “End”

I noticed that when I uncomment the first part (lines 7-22), I get the error:

outRad = GPF.createProduct(“Rad2Refl”, parRad, outResam)
RuntimeError: org.esa.snap.core.gpf.OperatorException: Sensor ‘SLSTR 500m’ not compliant with input product. Please check your selection.

while when they are commented out everything is working fine. I’m not sure it should work like that. I’m not changing the input, am I?

2. The above code is part of two scripts of mine one that is almost identical to the above and another one a little longer (it has a few other features as well). Yesterday evening I run both, the smaller one worked fine and the longer one gave me the error mentioned above. I tried several times both and I checked that the part of the 4 operators was identical between the scripts. This morning I run them again, without changing anything and none of them worked. I find this strange and it reminded me this topic: Best practice to convert and reproject Sentinel-3 radiances to reflectance, talking, unfortunately, about the same operator (Rad2Ref), although regarding a different aspect. Marpet is telling somewhere that there was a bug presented randomly and I thought it could apply to my case too. Just in case, I run snappy on windows and I have recently installed it so I have SNAP version 6.

3. I’ve noticed that Subset and Reproject don’t work directly with the .xml (which contains bands of different spatial resolution), without resampling first. It would be more practical if they worked even on mixed products so that people wouldn’t have to repeat the operators twice, once for the 500m and once for the 1Km bands…

4. I’ve tried two different workflows on snappy:
   Resample -> Rad2Ref -> Subset -> Reprojet -> Subset

and

Resample -> Reproject ->Rad2Ref -> Subset

and although none gives me an error, the second one gives me wrong results. The radiances haven’t been divided by 10.000 and even if I manually divide the result, some points still present differences that go from -0.7077 to 0.6235 (quite high for a variable that should stay within [0,1] (reflectance)).

Finished!! I’m sorry for the giant post. I just wanted to share my probably insignificant observations with everyone, in case snappy has made others too getting doubts on their mental health…

Have a nice day.

[aksshre]

can you just replace your inp with this
band_type = ‘Sen3_SLSTRL1B_500m’
reader1000 = snappy.ProductIO.getProductReader(band_type)
in_file = snappy.File("/MyPath/S3B_SL_1_RBT____20190219T092217_20190219T092517_20190219T113119_0179_022_150_2160_MAR_O_NR_003.SEN3/xfdumanifest.xml")
inp = reader1000.readProductNodes(in_file, None)

By default readProduct will read Sen3_SLSTRL1B_1km. This is 1 Km product and will only have band information for S7-S9( and F1,F2).

Also i dnot know why you are doing the resampling, by default, for 500m product, the width and height is 3000,2400.

Let me know if it works.
Thanks
Ashish

[Nobody]

Thank you for your answer. I’ll try it and let you know.

I don’t know why you say it reads only tyhe 1Km bands… I guess you refer to my first post but I 'd like to inform you that later it saw all the bands (286) without me changing anything…
Currently I work on the script posted, which is loading all the bands (from s1 to s9and the f ones). The reason I’m using resampling is that neither Rad2Ref nor Reproject work with the mixed band resolution product I’m reading. I know it seems stupid to resample the S1-S6 bands too (which are already at 500m) and it is stupid transforming the thermal bands into radiances (just because after subseting neither Rad2Ref nor Reproject works), but it’s the only way I managed to make things work. Of course for the thermal bands I will have to run a completely new procedure from the beginning, since the output of the above (for the thermal bands) is completely wrong and useless.
However, I will try your procedure and check whether it manages to load only the 500m bands.

[Nobody]

First I’d like to correct my previous post. Rad2Refl doesn’t nedd resampling, since it works fine on a mixed resolution product (with both 500m and 1km bands). The operators that give me error are Reproject and Subset.

I just tried the code you suggested and it gives me 167 bands, among which S1-S9 and the f ones. It looks like all the bands are loaded at 500m (the thermal too), so there is no need for resampling anymore, and the operators Reproject and Subset work perfectly. The Rad2Refl though still gives me the error: RuntimeError:
org.esa.snap.core.gpf.OperatorException: Sensor ‘SLSTR 500m’ not compliant with input product. Please check your selection.
when I call it after both Reproject and Subset for the s1 to s6 an bands, and still gives me wrong results without error after only the Reproject operator (see point 4 of my second post).

[marpet]

It seems you have worked it out, almost.
There is an issue with using the SLSTR data from the command line.
I’ve noted it here:
https://senbox.atlassian.net/browse/SIIITBX-238

You should explicitly define which reader you want, as there three available. Ashish mentioned this already.

That Rad2Refl is complaining about the input product, might be caused by a SLSTR product in 1KM resolution.

[aksshre]

If i understand your program correctly, you just need reflectance for band S1-S6 (an). As i and marpet told, you should explicitly define reader, in your case it should be Sen3_SLSTRL1B_500m. So i copied your program and try to run with just imp replaced. I did not get any error.

import snappy
from snappy import GPF, jpy

band_type = ‘Sen3_SLSTRL1B_500m’
reader1000 = snappy.ProductIO.getProductReader(band_type)
in_file = snappy.File("/file/Sentinel3-a/S3A_SL_1_RBT____20170210T130027_20170210T130327_20170211T211143_0179_014_152_3240_LN2_O_NT_002.SEN3/xfdumanifest.xml")
inp = reader1000.readProductNodes(in_file, None)

#resample
HashMap = jpy.get_type(‘java.util.HashMap’)
parResam = HashMap()
parResam.put(“targetWidth”, 3000)
parResam.put(“targetHeight”, 2400)
parResam.put(“upsampling”, “Nearest”)
parResam.put(“downsampling”, “First”)
parResam.put(“flagDownsampling”, “First”)
parResam.put(“resampleOnPyramidLevels”, “true”)

#radtoref
HashMap = jpy.get_type(‘java.util.HashMap’)
parRad = HashMap()
parRad.put(“sensor”, “SLSTR_500m”)
parRad.put(“conversionMode”, “RAD_TO_REFL”)
parRad.put(“copyTiePointGrids”, “false”)
parRad.put(“copyFlagBandsAndMasks”, “false”)
parRad.put(“copyNonSpectralBands”, “false”)

#reproject
parProj = HashMap()
parProj.put(“crs”, ‘PROJCS[“WGS 84 / UTM zone 32N”, GEOGCS[“WGS 84”, DATUM[“World Geodetic System 1984”, SPHEROID[“WGS 84”, 6378137.0, 298.257223563, AUTHORITY[“EPSG”,“7030”]], AUTHORITY[“EPSG”,“6326”]], PRIMEM[“Greenwich”, 0.0, AUTHORITY[“EPSG”,“8901”]], UNIT[“degree”, 0.017453292519943295], AXIS[“Geodetic longitude”, EAST], AXIS[“Geodetic latitude”, NORTH], AUTHORITY[“EPSG”,“4326”]], PROJECTION[“Transverse_Mercator”, AUTHORITY[“EPSG”,“9807”]], PARAMETER[“central_meridian”, 9.0], PARAMETER[“latitude_of_origin”, 0.0], PARAMETER[“scale_factor”, 0.9996], PARAMETER[“false_easting”, 500000.0], PARAMETER[“false_northing”, 0.0], UNIT[“m”, 1.0], AXIS[“Easting”, EAST], AXIS[“Northing”, NORTH], AUTHORITY[“EPSG”,“32632”]]’)
parProj.put(“resampling”, “Nearest”)
parProj.put(“orthorectify”, “false”)
parProj.put(“noDataValue”, -32768.0)
parProj.put(“includeTiePointGrids”, “true”)
parProj.put(“addDeltaBands”, “false”)

#subset
parSub1 = HashMap()
parSub1.put(“bandNames”, “S1_reflectance_an,S2_reflectance_an,S3_reflectance_an,S5_reflectance_an,S6_reflectance_an”)
parSub1.put(“copyMetadata”, “true”)

outResam = GPF.createProduct(“Resample”, parResam, inp)
outRad = GPF.createProduct(“Rad2Refl”, parRad, outResam)
outSub1 = GPF.createProduct(“Subset”, parSub1, outRad)
outProj = GPF.createProduct(“Reproject”, parProj, outSub1)
outSub2 = GPF.createProduct(“Subset”, parSub1, outProj)

[Nobody]

Ok, maybe I wasn’t clear enough.

My second post (the one with the code) was already working (when lines 7-22 were commented out).
My post was just about some considerations of mine regarding the correct functioning of snappy (which probably works as it should, but it wasn’t so clear to me so I preferred writing it in case someone else had my doubts too).
Ashish said that not specifying the resolution would by default load only the 1km bands (or maybe all the bands at 1km) which is not what I had noticed, and suggested me to clearly specify the resolution I wanted to load, which I did and it worked fine and I posted my results as asked.

That’s all. I don’t face any problem running my code and have never used GPT.
Sorry for the confusion and thanks for your time.