33import com .bc .ceres .binding .ValueRange ;
44import com .bc .ceres .core .ProgressMonitor ;
55import eu .esa .opt .processor .rad2refl .Rad2ReflConstants ;
6- import org .esa .snap .core .datamodel .Band ;
7- import org .esa .snap .core .datamodel .FlagCoding ;
8- import org .esa .snap .core .datamodel .GeoCoding ;
9- import org .esa .snap .core .datamodel .GeoPos ;
10- import org .esa .snap .core .datamodel .Product ;
11- import org .esa .snap .core .datamodel .ProductData ;
6+ import org .esa .snap .core .datamodel .*;
127import org .esa .snap .core .gpf .Operator ;
138import org .esa .snap .core .gpf .OperatorException ;
149import org .esa .snap .core .gpf .OperatorSpi ;
1813import org .esa .snap .core .gpf .annotations .SourceProduct ;
1914import org .esa .snap .core .gpf .annotations .TargetProduct ;
2015import org .esa .snap .core .util .ProductUtils ;
16+ import org .esa .snap .core .util .math .RsMathUtils ;
2117import org .esa .snap .idepix .core .IdepixConstants ;
2218import org .esa .snap .idepix .core .seaice .LakeSeaIceAuxdata ;
2319import org .esa .snap .idepix .core .seaice .LakeSeaIceClassification ;
@@ -104,6 +100,11 @@ public class IdepixOlciClassificationOp extends Operator {
104100 "Slower, but in general more precise." )
105101 private boolean useSrtmLandWaterMask ;
106102
103+ @ Parameter (defaultValue = "true" ,
104+ label = " If selected, harmonized radiances are taken from O2 product and used in NN input" ,
105+ description = " If selected, harmonized radiances are taken from O2 product and used in NN input " )
106+ private boolean useO2HarmonizedRadiancesForNN ;
107+
107108
108109 @ SourceProduct (alias = "l1b" , description = "The L1b product." )
109110 private Product l1bProduct ;
@@ -124,13 +125,23 @@ public class IdepixOlciClassificationOp extends Operator {
124125
125126 private Band [] olciReflBands ;
126127
128+ private Band solarFlux13Band ;
129+ private Band solarFlux14Band ;
130+ private Band solarFlux15Band ;
131+
132+ private RasterDataNode szaBand ;
133+
134+ private Band harmoRad13Band ;
135+ private Band harmoRad14Band ;
136+ private Band harmoRad15Band ;
137+
127138 private Band surface13Band ;
128139 private Band trans13Band ;
129140
130141 private static final String OLCI_2018_NET_NAME = "11x10x4x3x2_207.9.net" ;
131142 private static final String OLCI_2018_NN_THRESHOLDS_FILE = "11x10x4x3x2_207.9-thresholds.json" ;
132- private static final String OLCI_202306_NET_NAME = "class-sequential-i21x42x8x4x2o1-5489.net" ;
133- private static final String OLCI_202306_NN_THRESHOLDS_FILE = "class-sequential-i21x42x8x4x2o1-5489-thresholds.json" ;
143+ private String OLCI_202306_NET_NAME = "class-sequential-i21x42x8x4x2o1-5489.net" ;
144+ private String OLCI_202306_NN_THRESHOLDS_FILE = "class-sequential-i21x42x8x4x2o1-5489-thresholds.json" ;
134145
135146 private static final double THRESH_LAND_MINBRIGHT1 = 0.3 ;
136147 private static final double THRESH_LAND_MINBRIGHT2 = 0.25 ; // test OD 20170411
@@ -154,6 +165,12 @@ public class IdepixOlciClassificationOp extends Operator {
154165 @ Override
155166 public void initialize () throws OperatorException {
156167 setBands ();
168+
169+ OLCI_202306_NET_NAME = useO2HarmonizedRadiancesForNN ? "class-sequential-i21x42x8x4x2o1-5489-new-lake-ice-o2harm.net" :
170+ "class-sequential-i21x42x8x4x2o1-5489.net" ;
171+ OLCI_202306_NN_THRESHOLDS_FILE =
172+ useO2HarmonizedRadiancesForNN ? "class-sequential-i21x42x8x4x2o1-5489-new-lake-ice-o2harm-thresholds.json" :
173+ "class-sequential-i21x42x8x4x2o1-5489-thresholds.json" ;
157174 readSchillerNeuralNets ();
158175 readNNThresholds ();
159176 nnInterpreter = IdepixOlciCloudNNInterpreter .create ();
@@ -171,6 +188,15 @@ public void initialize() throws OperatorException {
171188 initLakeSeaIceClassification ();
172189
173190 if (o2CorrProduct != null ) {
191+ szaBand = l1bProduct .getRasterDataNode ("SZA" );
192+
193+ solarFlux13Band = l1bProduct .getBand ("solar_flux_band_13" );
194+ solarFlux14Band = l1bProduct .getBand ("solar_flux_band_14" );
195+ solarFlux15Band = l1bProduct .getBand ("solar_flux_band_15" );
196+
197+ harmoRad13Band = l1bProduct .getBand ("radiance_13" );
198+ harmoRad14Band = l1bProduct .getBand ("radiance_14" );
199+ harmoRad15Band = l1bProduct .getBand ("radiance_15" );
174200 surface13Band = o2CorrProduct .getBand ("surface_13" );
175201 trans13Band = o2CorrProduct .getBand ("trans_13" );
176202 gf = new GeometryFactory ();
@@ -202,11 +228,11 @@ void readNNThresholds() {
202228 for (NNThreshold t : NNThreshold .values ()) {
203229 if (m .containsKey (t .name ())) {
204230 t .range = new ValueRange ((Double ) ((JSONArray ) m .get (t .name ())).get (0 ),
205- (Double ) ((JSONArray ) m .get (t .name ())).get (1 ),
206- true ,
207- false );
231+ (Double ) ((JSONArray ) m .get (t .name ())).get (1 ),
232+ true ,
233+ false );
208234 } else {
209- t .range = new ValueRange (0.0 , 0.0 ,true , false );
235+ t .range = new ValueRange (0.0 , 0.0 , true , false );
210236 }
211237 }
212238 } catch (FileNotFoundException e ) {
@@ -274,6 +300,39 @@ public void computeTileStack(Map<Band, Tile> targetTiles, Rectangle rectangle, P
274300 trans13Tile = getSourceTile (trans13Band , rectangle );
275301 }
276302
303+ Tile szaTile = null ;
304+ if (szaBand != null ) {
305+ szaTile = getSourceTile (szaBand , rectangle );
306+ }
307+
308+ Tile solarFlux13Tile = null ;
309+ if (solarFlux13Band != null ) {
310+ solarFlux13Tile = getSourceTile (solarFlux13Band , rectangle );
311+ }
312+ Tile solarFlux14Tile = null ;
313+ if (solarFlux14Band != null ) {
314+ solarFlux14Tile = getSourceTile (solarFlux14Band , rectangle );
315+ }
316+ Tile solarFlux15Tile = null ;
317+ if (solarFlux15Band != null ) {
318+ solarFlux15Tile = getSourceTile (solarFlux15Band , rectangle );
319+ }
320+
321+ Tile harmoRad13Tile = null ;
322+ if (harmoRad13Band != null ) {
323+ harmoRad13Tile = getSourceTile (harmoRad13Band , rectangle );
324+ }
325+
326+ Tile harmoRad14Tile = null ;
327+ if (harmoRad14Band != null ) {
328+ harmoRad14Tile = getSourceTile (harmoRad14Band , rectangle );
329+ }
330+
331+ Tile harmoRad15Tile = null ;
332+ if (harmoRad15Band != null ) {
333+ harmoRad15Tile = getSourceTile (harmoRad15Band , rectangle );
334+ }
335+
277336 final Band olciQualityFlagBand = l1bProduct .getBand (IdepixOlciConstants .OLCI_QUALITY_FLAGS_BAND_NAME );
278337 final Tile olciQualityFlagTile = getSourceTile (olciQualityFlagBand , rectangle );
279338
@@ -315,10 +374,16 @@ public void computeTileStack(Map<Band, Tile> targetTiles, Rectangle rectangle, P
315374 // todo: for cglops, coastlines are treated as LAND
316375 if ((isLandFromAppliedMask && !isInlandWaterFromAppliedMask ) || isCoastlineFromAppliedMask ) {
317376 classifyOverLand (olciReflectanceTiles , cloudFlagTargetTile , nnTargetTile ,
318- surface13Tile , trans13Tile , x , y );
377+ surface13Tile , trans13Tile ,
378+ solarFlux13Tile , solarFlux14Tile , solarFlux15Tile ,
379+ harmoRad13Tile , harmoRad14Tile , harmoRad15Tile ,
380+ szaTile , x , y );
319381 } else {
320382 classifyOverWater (olciQualityFlagTile , olciReflectanceTiles ,
321- cloudFlagTargetTile , nnTargetTile , x , y , isInlandWaterFromAppliedMask );
383+ cloudFlagTargetTile , nnTargetTile ,
384+ solarFlux13Tile , solarFlux14Tile , solarFlux15Tile ,
385+ harmoRad13Tile , harmoRad14Tile , harmoRad15Tile ,
386+ szaTile , x , y , isInlandWaterFromAppliedMask );
322387 }
323388 }
324389 }
@@ -334,10 +399,16 @@ private boolean classifyCoastline(Tile olciQualityFlagTile, int x, int y, int wa
334399 }
335400
336401 private void classifyOverWater (Tile olciQualityFlagTile , Tile [] olciReflectanceTiles ,
337- Tile cloudFlagTargetTile , Tile nnTargetTile , int x , int y , boolean isInlandWater ) {
402+ Tile cloudFlagTargetTile , Tile nnTargetTile ,
403+ Tile solarFlux13Tile , Tile solarFlux14Tile , Tile solarFlux15Tile ,
404+ Tile harmoRad13Tile , Tile harmoRad14Tile , Tile harmoRad15Tile ,
405+ Tile szaTile , int x , int y , boolean isInlandWater ) {
406+
407+ final double [] o2HarmoRefls = getO2HarmoReflectances (solarFlux13Tile , solarFlux14Tile , solarFlux15Tile ,
408+ harmoRad13Tile , harmoRad14Tile , harmoRad15Tile , szaTile , x , y );
338409
410+ final double nnOutput = getOlciNNOutput (x , y , olciReflectanceTiles , o2HarmoRefls );
339411
340- double nnOutput = getOlciNNOutput (x , y , olciReflectanceTiles );
341412 if (!cloudFlagTargetTile .getSampleBit (x , y , IdepixConstants .IDEPIX_INVALID )) {
342413 cloudFlagTargetTile .setSample (x , y , IdepixConstants .IDEPIX_CLOUD_AMBIGUOUS , false );
343414 cloudFlagTargetTile .setSample (x , y , IdepixConstants .IDEPIX_CLOUD_SURE , false );
@@ -387,9 +458,14 @@ private void classifyOverWater(Tile olciQualityFlagTile, Tile[] olciReflectanceT
387458 private void classifyOverLand (Tile [] olciReflectanceTiles ,
388459 Tile cloudFlagTargetTile , Tile nnTargetTile ,
389460 Tile surface13Tile , Tile trans13Tile ,
390- int x , int y ) {
461+ Tile solarFlux13Tile , Tile solarFlux14Tile , Tile solarFlux15Tile ,
462+ Tile harmoRad13Tile , Tile harmoRad14Tile , Tile harmoRad15Tile ,
463+ Tile szaTile , int x , int y ) {
391464
392- final double nnOutput = getOlciNNOutput (x , y , olciReflectanceTiles );
465+ final double [] o2HarmoRefls = getO2HarmoReflectances (solarFlux13Tile , solarFlux14Tile , solarFlux15Tile ,
466+ harmoRad13Tile , harmoRad14Tile , harmoRad15Tile , szaTile , x , y );
467+
468+ final double nnOutput = getOlciNNOutput (x , y , olciReflectanceTiles , o2HarmoRefls );
393469
394470 if (!cloudFlagTargetTile .getSampleBit (x , y , IdepixConstants .IDEPIX_INVALID )) {
395471 cloudFlagTargetTile .setSample (x , y , IdepixConstants .IDEPIX_CLOUD_AMBIGUOUS , false );
@@ -448,6 +524,32 @@ private void classifyOverLand(Tile[] olciReflectanceTiles,
448524 }
449525 }
450526
527+ private static double [] getO2HarmoReflectances (Tile solarFlux13Tile , Tile solarFlux14Tile , Tile solarFlux15Tile , Tile harmoRad13Tile , Tile harmoRad14Tile , Tile harmoRad15Tile , Tile szaTile , int x , int y ) {
528+ float sza ;
529+ float sf13 ;
530+ float sf14 ;
531+ float sf15 ;
532+ float harmoRad13 ;
533+ float harmoRad14 ;
534+ float harmoRad15 ;
535+ if (szaTile != null && solarFlux13Tile != null && solarFlux14Tile != null && solarFlux15Tile != null &&
536+ harmoRad13Tile != null && harmoRad14Tile != null && harmoRad15Tile != null ) {
537+ sza = szaTile .getSampleFloat (x , y );
538+ sf13 = solarFlux13Tile .getSampleFloat (x , y );
539+ sf14 = solarFlux14Tile .getSampleFloat (x , y );
540+ sf15 = solarFlux15Tile .getSampleFloat (x , y );
541+ harmoRad13 = harmoRad13Tile .getSampleFloat (x , y );
542+ harmoRad14 = harmoRad14Tile .getSampleFloat (x , y );
543+ harmoRad15 = harmoRad15Tile .getSampleFloat (x , y );
544+
545+ final double harmoRefl13 = RsMathUtils .radianceToReflectance (harmoRad13 , sza , sf13 );
546+ final double harmoRefl14 = RsMathUtils .radianceToReflectance (harmoRad14 , sza , sf14 );
547+ final double harmoRefl15 = RsMathUtils .radianceToReflectance (harmoRad15 , sza , sf15 );
548+ return new double []{harmoRefl13 , harmoRefl14 , harmoRefl15 };
549+ }
550+ return null ;
551+ }
552+
451553 private boolean isOlciLandPixel (int x , int y , Tile olciL1bFlagTile , int waterFraction ) {
452554 if (waterFraction < 0 ) {
453555 boolean landFlag = olciL1bFlagTile .getSampleBit (x , y , IdepixOlciConstants .L1_F_LAND );
@@ -495,7 +597,7 @@ private boolean isOlciInlandWaterPixel(int x, int y, Tile olciL1bFlagTile, int w
495597 }
496598 }
497599
498- private double getOlciNNOutput (int x , int y , Tile [] rhoToaTiles ) {
600+ private double getOlciNNOutput (int x , int y , Tile [] rhoToaTiles , double [] o2HarmoRefls ) {
499601 SchillerNeuralNetWrapper nnWrapper ;
500602 try {
501603 nnWrapper = olciAllNeuralNet .get ();
@@ -506,6 +608,11 @@ private double getOlciNNOutput(int x, int y, Tile[] rhoToaTiles) {
506608 for (int i = 0 ; i < nnInput .length ; i ++) {
507609 nnInput [i ] = Math .sqrt (rhoToaTiles [i ].getSampleFloat (x , y ));
508610 }
611+ if (o2HarmoRefls != null ) {
612+ for (int i = 0 ; i < 3 ; i ++) {
613+ nnInput [i + 12 ] = Math .sqrt (o2HarmoRefls [i ]);
614+ }
615+ }
509616 return nnWrapper .getNeuralNet ().calc (nnInput )[0 ];
510617 }
511618
0 commit comments