Dear Nasir Hameed,

two things come to my mind:

1.) While Band 8 and Band 4 of Sentinel-2 have a spatial resolution of 10m, Band 7 comes at 20m. This will affect the calculated values of NDVIre3 e.g. at the boundaries of different landcover types. I did not read the linked publication, but maybe it is also appropriate to use band 8A (which comes at 20m) instead of Band 8 to avoid differences in spatial resolution?

2.) Furthermore, Band 8 and Band 7 record at a more similar wavelength (- in fact they even overlap a bit -) than e.g. Band 8 and Band 4. Differences between the signal in Band 8 and Band 7 are therefore expected to be comparably low. Consequently, the remaining noise of each channel might have a larger impact on the calculation of NDVIre3 than on the the calculation of NDVI.

Hope that helps,

Sebastian

I am thankful to you for your input, please note that;

The SCP plugin of QGIS converts the Sentinel-2 images (Level-1C) to the physical measure of Top Of Atmosphere reflectance (TOA), applying a simple atmospheric correction using the DOS1 method (Dark Object Subtraction 1), which is an image-based technique. Bands with 20m spatial resolution are resampled to 10m resolution without changing the original pixel value (i.e. one 20m pixel is divided in four 10m pixels with the same value).

I am interested to compute NDVIre3 for which Band 8 and 7 of this sensor are involved. While, to calculate NDVIre3n (red-edge 3 narrow) Band 8a and 7 are used i.e. (B8a - B7) / (B8a + B7).

Unfortunately, both indices produce 'noisy' raster having discrepancies in output range of values and appearance clearly evident when compared to existing land cover through simple visual inspection.

Have you tried to download Sentinel 2 - level 2A from Copernicus data Hub and then calculate NDVIre3on SNAP software?

Marco Sozzi, ESA is currently providing Level-2A products for Europe only.

https://earth.esa.int/web/sentinel/user-guides/sentinel-2-msi/product-types/level-2a

The L1C products can be converted to Level-2A to retrieve scaled reflectance values using Sen2Cor external plugin processor in SNAP, . In addition to Sen2Cor and QGIS, I have derived surface reflectance values using multiple software but still problem persists.

Thank you Nasir Hameed, I didn't that Level-2A is not available outside Europe. However, the algorithms used by ESA and Sen2Cor are different.

On SNAP it should work. I tried to calculate NDVIre3 in the attached file.

Marco Sozzi I am thankful to you for providing me assistance in this regard.

I am trying to figure out possible source of problem in SNAP and revisiting pre-processing steps.

A glance on your raster reveals that higher values are concentrated along the agricultural field edges rather than middle of the fields if cultivated otherwise!

I have revisited pre-processing stage. Performed atmospheric corrections on L1C data using Sen2Cor in SNAP and retrieved L2A product. It didn't made any difference. The indices are again noisy with weird looking values.

I do not know if this is the case, but computing NDVI from atmospherically corrected products, I often experienced problematic output, when the bands contained pixels with negative values - these can occur as result of over-correction or noise and artifacts in the images in the dark regions. Using bands with negative values in any normalized difference formula then leads to occurrence of values way outside of the expected normalized difference values. Check your input bands for negative and zero reflectance values and the resulting NDVIre3 for values outside of .

Quick fix is just to turn any negative values to zero, or better (to avoid still possible division by zero) even to a very low positive reflectance value before computing the vegetation index. In case of quantized reflectance in L2A Sentinel-2, where reflectance 1.0 is represented by DN=10000, I would turn negative and zero values to DN=1 or reflectance 0.0001.

Of course in case of general over-correction with more than few scattered negative values in the image it is needed to look at the atmospheric correction and why it occurs.

And in case the image noise is the reason why it happens, an image smoothing technique would help, although sacrifying some of the detail level.

Tomáš Brunclík Thank you for your feedback. I had applied DOS1 atmospheric correction algorithm available in QGIS to ensure there are no pixels with zero data value and the pixels size is uniform for both input bands (B7, B8) prior to compute NDVIre3 & NDVIre3n. Subsequently, I have carefully analyzed the histogram of scaled output product to ensure there are no negative/zero values within the input bands. However, the output product contains unrealistic data range and associated noise. The resultant noise is not exclusively limited to hill slopes with shadows. Image filtering techniques like smoothing were found to be little helpful in this context.

Secondly, I derived scaled ground reflectance product L2A using Sen2Cor plugin in SNAP but encountered almost same results. The noise is again well distributed across the extent. Performing image smoothing techniques didn't produced realistic results thereby.

Thirdly, I experimented with different raster resampling approaches as advised by the STEP forum but no avail.

Could you please work on a subset of S2 raster to ascertain the source of this issue? I will be grateful for this support.

Nasir Hameed I tried to compute NDVIre3 on a subset of product

S2A_MSIL2A_20180906T101021_N0208_R022_T33UWR_20180906T131549.SAFE. It is Level-2 product distributed by ESA, but it should be similar to product obtained using SEN2COR correction.

• First let us look at the selected area in 11-6-2 color compostion (img1). It is from 20m composite image of bands 2 3 4 5 6 7 11 12 8A i created of the region of interest. The bands are those output by SEN2COR in 20m resolution.
• I then added individual JP2 files for bands 07 and 08, did spatial subset of these, redicing resolution of band 08 to 20m. I also set value 0 to NODATA value in the resulting map layers. From these bands I computed NDVIre3 using the formula mentioned in your first post in the thread. In the map it looks weird indeed (img2), there is little contrast between various surface classes, but it has seemingly enhanced in 3D the borders in some features - this is probably because of misalignment of the bands 7 and 8 which was apparent already before resampling band 8 to 20m. The range of values is approx -0.984 to 0.781.
• Then I tried another thing: I computed NDVIre3 from band 07 and 8A (8A instead of band 08) - this band is already 20m resolution - the same as 07. It is narrow band in the same part of spectrum as broad band 08. The result is in image (img3). This looks a little better - no 3D effect and there is distinction between various types of surfaces. The range is approx -0.915 to 0.610. I do not know if using this band instead of 08 is right, but generally, I would suppose it should be better for doing differencing with band 07 (normalized difference), because band 08 is so wide, that it even overlaps part of spectral region of band 07. Band 8A is spectrally well apart and approximately as wide as 07.
• For better comparison of the two NDVIre3 variants I stretched the gray pallete of the 8A based one, and copied the style to the 08 one. This way both have enhanced contrast and the same absolute value corresponds to the same grey level in both: (img4) - NDVIre3[8A], (img5) - NDVIre3[08]. Comparing these two, the [8A] better discriminates various land cover types, but is noisy. [08] is not so noisy, but is affected by misalignment at the edge of land cover types.
• I tried if smoothing can mitigate the noise in NDVIre3[8A] by smoothing the input bands 07 and 8A with GRASS tools in QGIS. I used the tool r.mfilter with a filter doing simple 3x3 pixels averaging, then used QGIS map calculator on the exported smoothed bands 07 and 8A. It works, although some noise remains around edges and of course, detail level decreases (img6).

Try to compare the steps you did with mine, if there are any differences. For example making sure there are no negative values is not sufficient by looking at the histogram, because you can miss few negative pixels in it, unless you can set logarithmic scale in the counts. In QGIS it is sufficient to unstretch the colors by stretching the colors/grays to the full scale of the dataset and you see the actual range in the list of map layers. If the 08 band must be used and you also have problems with alignment, realignment would be necessary prior to computing the index.

@ Tomáš Brunclík I am really grateful for your anticipation, effort and time.

I have analyzed your results and followed steps you have taken to derive these indices. A closer look on attached images reveals that you are facing same problem I am encountering with. Certainly, as you have rightly pointed out that we have to eliminate 'misalignment/noise' in these bands prior to compute these indices. However, source of this noise/misalignment is not yet clear. May be we have to revisit 'resampling' of these bands for consistency.

In both cases, i.e. NDVIre3 and NDVIre3n, the higher values should be concentrated over areas of healthy vegetation as seen in the true color composite (img1.png). The range of values is still unrealistic with respect to ground. E.g. higher values can be noticed over uncultivated fields in the produced results.

Secondly, the L2A product distributed by ESA contains 3D effect in processed bands surely because of DEM employed in its atmospheric correction model. It induces errors which are quite visible in the resultant images.

Nasir Hameed, by the 3D effect I meant the light and dark effect on feature edges, which is due to misalignment of input bands. From comparing the bands in detail it seems, that while resampling might have added to this effect, some slight misalignment is apparent if you compare bands with 10m vs. 20m resolution even before resampling to common resolution. It could very well be problem of the software withe data format, not of the Sentinel data itself. I tried to view and resample the same product to 20m and then compute the NDVIre3 from bands 8 and 7 in SNAP, the software designed to work with the Sentinels data specifically, and there is no 3D effect (except on clouds, but these are shifted between bands due to paralax normally). The misalignment in QGIS and possibly other software may be because of the fact, that different geographic software handles differently the spatial registration metadata within the raster product or band. For example if a raster is registered into a Cartesian spatial reference system by corner pixels coordinates, it may be meant, that these coordinates are either centers of the corner pixels, or it can mean outward corners of the corner pixels. If the software handles it the other way around than it was meant in the band metadata, you would normally hardly notice a difference, except if you work with bands with different resolution - then a slight misalignment may show. But the noise is still there and also, the vegetation does not stand out from built-up areas (NDVIre3_0708_SNAP.png). I did not have the time to read about the NDVIre3, but perhaps this vegetation index may not be intended to delineate vegetation (as does NDVI) but to show certain differences within the vegetation cover?

Anyway, to mitigate the problem of misalignment, SNAP seems to be the right tool. Also filtering the bands before computing the index is possible (NDVIre3F_0708_SNAP.png).

Tomáš Brunclík Regarding 3D effect, I was referring to img2 and img5.png.

These indices have been used in many applications e.g.

Article SENTINEL-2A red-edge spectral indices suitability for discri...

https://www.sciencedirect.com/science/article/abs/pii/S0924271616306128?via%3Dihub

I had recently downloaded L2A product being distributed by ESA

S2A_MSIL2A_20181227T055231_N0211_R048_T43SCU_20181227T075542

and found discrepancies induced by atmospheric correction employed to get L2A product. The 3D effect is quite visible in resultant individual bands (attached).

However, I will revisit to ensure misalignment issues are eliminated prior to compute NDVIre3 and NDVIre3n.

I am indebted to your continuous guidance and help in this regard, please.

Nasir Hameed the hill shading effect seen in your images should be removed or strongly suppressed by the DEM employed during atmospheric correction, not caused by it. However, when you do the correction yourself with Sen2Cor, you have to make sure the terrain correction with DEM is used, because the default is to NOT use any DEM. You have to edit the GIPP config file and make sure the DEM is successfully downloaded the first time it is used. See chapter 2.7.7 in

http://step.esa.int/thirdparties/sen2cor/2.5.5/docs/S2-PDGS-MPC-L2A-SUM-V2.5.5_V2.pdf

Also some level-2 products distributed by ESA may be produced without employing terrain correction, especially older products. You can see if and which DEM was used in the metadata file MTD_MSIL2A.xml in the .SAFE directory of the product, just search for "dem" in the file (the key is PRODUCTION_DEM_TYPE). If there is a name of dem used (in my case DTED), the product should be terrain corrected.

My area of interest is quite flat, so I do not know how well should the hill shading be suppressed, but looks to me like it is not corrected at all in your images. It should be less visible, than in the corresponding level-1 product stretched to the same contrast.

Tomáš Brunclík Thank you for your feedback and help in this context. I am working to resolve all of these imperfections prior to get reliable results.