Hello,
I am trying to extract pixel values (Raster->Export->Extract Pixel Values) from many different Sentinel OLCI .xml files (many source paths). I add my coordinates and I select my output directory and then I select extract. Then a successful message appears (The pixel extraction tool has run successfully and written the result files to C:…etc), but there is no an output file (.txt *file /prefix: pixEx ) in the output directory.
Can you help me?
Thank you a lot,
Anastasia
This often happens because some constraints are not met. Either the allowed time distance too high or the coordinates are not located in the source products or the source products are rejected for some other reason.
Unfortunately the error reporting within the desktop application is not optimal. When you invoke it from the command line you should see the reasons why no pixels are extracted.
I’ve attached a template file where you can insert the parameters from the GUI.
PixExTemplate.xml (409 Bytes)
On the command line you can invoke it like this:
gpt PixExTemplate.xml
Hello there @marpet, first of all thanks for making yourself available and dedicating your time helping us all over the forum.
I am having a similar issue, so I thought it would be ok to post in the same thread.
I have a region of interest in ESRI Shapefile format, I can subset my image using it, ok. But how do I extract all the values from the pixels within my shapefile? I don’t see any options related to it in the GUI…
I have managed to programmatically extract a single pixel by a given coordinate using python, but I am having a hard time manupulating the Sentinel-3 NetCDF files. Any directions would be of a great help on this subject.
Hi, and thanks for the nice words.
The question is do you want the extracted data in an ASCII format or or do you still want an image format? It is also different if you want to do it on a batch of products or just a few in the desktop.
And also thanks for the quick reply, you see? You deserve more than nice words for helping us with such complex tasks.
Both will work, so long I can get rid of all the data outside the shapefile/ROI or get all the data inside it.
My final goal is to have the value of every individual pixel within the green area for all available bands of a given Sentinel-3 WFR or EFR product, even if that data is invalid or under a cloud pixel.
My biggest struggle so far is that the individual Oa**_reflectance.nc files inside the image folder seems not to have geographical information attached to them, and that information is inside the geo_coordinates.nc
I have tried converting the WFR NetCDF images into a GeoTIFF/BigTIFF to make them easier to manipulate but the resulting stack was like 20Gb+ and I needed only to operate in the pixels within this green shapefile above…
(just for the sake of curiosity, this is part of the Manacapuru river in Brazil btw)
When you say ‘different’ do you mean the method? I’ve managed to do it on 209 WFR images on a desktop, but I was extracting all the band values for a single pixel.
You are correct, I don’t really know how computationally expensive it will be when I do it for all the pixels within the shapefile instead of only one pixel…
This is my python3 code so far, but it is not very efficient because it loads up in memory all the NetCDF data from all the S3 bands of a given image:
bands_dictionary is a list of numpy N-D-Arrays containing the NetCDF data for ALL the bands
lon and lat are the longitude and latitude dimensions inside de file geo_coordinates.nc
target_lon and target_lat are the point of interest to be extracted like: -60.014493, -3.158980
def get_point_data_in_bands(self, bands_dictionary, lon=None, lat=None, target_lon=None, target_lat=None):
lat = lat[:, :, np.newaxis]
lon = lon[:, :, np.newaxis]
grid = np.concatenate([lat, lon], axis=2)
vector = np.array([target_lat, target_lon]).reshape(1, 1, -1)
subtraction = vector - grid
dist = np.linalg.norm(subtraction, axis=2)
result = np.where(dist == dist.min())
target_x_y = result[0][0], result[1][0]
rad_in_bands = []
for b in bands_dictionary:
rad = bands_dictionary[b][target_x_y]
rad_in_bands.append(rad)
return target_x_y, rad_in_bands
The output target_x_y is just in case I want to re-use it later, and represents a tuple containing the absolute X,Y position ex: (3825, 726) in the image matrix (like SNAP in the Pin Manager tool):
And rad_in_bands is what I actually want, a list with the radiance of the pixel in each band e.g.:
[-0.01853694,
-0.013849303,
0.004828036,
0.029529706,
0.043244734,
0.083712265,
0.11280863,
0.12397839,
0.12610249,
0.12782373,
0.11899777,
0.070143744,
0.07380597,
0.04839015,
0.042658776,
0.02198553]
I feel like I just did way too much code gymnastics to do something that probably has a better way of doing. And before I started writing more spaghetti code to do the same above for all the pixels inside a given ESRI shapefile:
I thought it would be smarter to ask you @marpet how to do it the correct way 
Ok, sorry for flooding the thread but I have just found this answer from you:
And it is exactly what I need, now, how do I do it from the command line?
Good that you have found this way. That’s what I wanted to explain for the Desktop use case.
When you use python you can use snappy. This helps reading the data.
https://senbox.atlassian.net/wiki/spaces/SNAP/pages/19300362/How+to+use+the+SNAP+API+from+Python
However, the simplest way to do on a batch of products is probably to use gpt. There are two wiki pages available which explain how to do it. They are linked here:
https://senbox.atlassian.net/wiki/spaces/SNAP/pages/70503053/Processing
For your use case, you need 3 operators chained in a graph.
I’ve attached a template. shapefileExtraction.xml (1.5 KB)
In the first step, the data is tailored to the bounding box. Here I just added a random rectangle.
In the second step, the shapefile is imported and in the last step, the data outside of the shapefile is masked.
Unfortunately, this approach does not export only the pixels covered by your shapefile but always the bounding box. Especially when exporting to CSV/ASCII format the amount of data is huge.
Thanks so much for all the tutorials and also for the xml template @marpet
through some trial and error I’ve managed to get snappy up and running system-wide alongside GDAL:
So…
Is there a way of making the gpt/gpf or snappy to simulate this behavior of the GUI?:
I did it in SNAP Desktop GUI and got this file: S3A_OL_2_WFR____20190309T141223_20190309T141523_20190310T211622_0179_042_167_3060_MAR_O_NT_002_manacapuru_Mask.txt (3.1 MB)
And I didn’t tested it to see where the lon/lat pixel data will fall on a map, was planning on doing it later today, but… My real question is: does that mean that the pixels that I extracted this way in the above file are actually everything inside the bounding box rectangle of the shapefile or only the ones actually inside the shapefile?
Sorry to keep bothering you guys with this, you already did such a great job by making all these great tools and data freely available to us.
Thanks again and keep safe! 
David
gotta love this SNAP man… The points extracted by the GUI fall right inside the ROI/SHP file:
Unfortunately, this is not possible to do it on a stack of data.
Only what I already wrote regarding gpt is possible.
We have this improvement already on the todo list (SNAP-828).
What you could do and if it is applicable for your use case, you could create a list of coordinates and use the Pixel Extraction tool (usable on a stack of products even in the GUI).
This way you wouldn’t get all pixels covered by your region, but maybe enough?
The data exported with the mask pixel extraction could be a starting point for this list.
Hello again @marpet. Thanks a lot for all the hints so far.
I have been writing some python code to exploit Sentinel-3 NetCDF and create time-series from selected areas using the GPT tool. Everything was working ok for a while and I’ve got some nice spectral curves.
Recently I have tried to crop some more images and then this happened:
15/09/2020 17:01:44 - Calling SNAP-GPT using command string:
"C:\Program Files\snap\bin\gpt.exe" "D:\S3\L2_WFR_subset\MANACAPURU\AUX_S3A_OL_2_WFR____20191104T135002_20191104T135302_20191124T134508_0179_051_124_3060_MAR_O_NT_002\gpt_manacapuru.xml" -f CSV -t "D:\S3\L2_WFR_subset\MANACAPURU\S3A_OL_2_WFR____20191104T135002_20191104T135302_20191124T134508_0179_051_124_3060_MAR_O_NT_002_subset.txt" -Ssource=D:\S3\L2_WFR\S3A_OL_2_WFR____20191104T135002_20191104T135302_20191124T134508_0179_051_124_3060_MAR_O_NT_002.SEN3
INFO: org.esa.snap.core.gpf.operators.tooladapter.ToolAdapterIO: Initializing external tool adapters
INFO: org.esa.snap.core.util.EngineVersionCheckActivator: Please check regularly for new updates for the best SNAP experience.
Executing processing graph
INFO: org.hsqldb.persist.Logger: dataFileCache open start
45%90% done.
org.esa.snap.core.gpf.OperatorException: Not able to write product file: 'D:\S3\L2_WFR_subset\MANACAPURU\S3A_OL_2_WFR____20191104T135002_20191104T135302_20191124T134508_0179_051_124_3060_MAR_O_NT_002_subset.txt'
at org.esa.snap.core.gpf.graph.GraphProcessor$GPFImagingListener.errorOccurred(GraphProcessor.java:363)
at com.sun.media.jai.util.SunTileScheduler.sendExceptionToListener(SunTileScheduler.java:1646)
at com.sun.media.jai.util.SunTileScheduler.scheduleTile(SunTileScheduler.java:921)
at javax.media.jai.OpImage.getTile(OpImage.java:1129)
at com.sun.media.jai.util.RequestJob.compute(SunTileScheduler.java:247)
at com.sun.media.jai.util.WorkerThread.run(SunTileScheduler.java:468)
Caused by: org.esa.snap.core.gpf.OperatorException: Not able to write product file: 'D:\S3\L2_WFR_subset\MANACAPURU\S3A_OL_2_WFR____20191104T135002_20191104T135302_20191124T134508_0179_051_124_3060_MAR_O_NT_002_subset.txt'
at org.esa.snap.core.gpf.common.WriteOp.computeTile(WriteOp.java:360)
at org.esa.snap.core.gpf.internal.OperatorImage.computeRect(OperatorImage.java:80)
at javax.media.jai.SourcelessOpImage.computeTile(SourcelessOpImage.java:137)
at com.sun.media.jai.util.SunTileScheduler.scheduleTile(SunTileScheduler.java:904)
... 3 more
Caused by: java.lang.ArrayStoreException
at java.lang.System.arraycopy(Native Method)
at org.esa.snap.csv.dataio.writer.CsvProductWriter.writeBandRasterData(CsvProductWriter.java:142)
at org.esa.snap.core.gpf.common.WriteOp.writeTileRow(WriteOp.java:408)
at org.esa.snap.core.gpf.common.WriteOp.computeTile(WriteOp.java:326)
... 6 more
SEVERE: org.esa.snap.core.util.SystemUtils$SnapImagingListener: JAI error occurred: 'Problem occurs when computing a tile by the owner.' at com.sun.media.jai.util.SunTileScheduler@e17dd41
org.esa.snap.core.gpf.OperatorException: Not able to write product file: 'D:\S3\L2_WFR_subset\MANACAPURU\S3A_OL_2_WFR____20191104T135002_20191104T135302_20191124T134508_0179_051_124_3060_MAR_O_NT_002_subset.txt'
at org.esa.snap.core.gpf.common.WriteOp.computeTile(WriteOp.java:360)
at org.esa.snap.core.gpf.internal.OperatorImage.computeRect(OperatorImage.java:80)
at javax.media.jai.SourcelessOpImage.computeTile(SourcelessOpImage.java:137)
at com.sun.media.jai.util.SunTileScheduler.scheduleTile(SunTileScheduler.java:904)
at javax.media.jai.OpImage.getTile(OpImage.java:1129)
at com.sun.media.jai.util.RequestJob.compute(SunTileScheduler.java:247)
at com.sun.media.jai.util.WorkerThread.run(SunTileScheduler.java:468)
Caused by: java.lang.ArrayStoreException
at java.lang.System.arraycopy(Native Method)
at org.esa.snap.csv.dataio.writer.CsvProductWriter.writeBandRasterData(CsvProductWriter.java:142)
at org.esa.snap.core.gpf.common.WriteOp.writeTileRow(WriteOp.java:408)
at org.esa.snap.core.gpf.common.WriteOp.computeTile(WriteOp.java:326)
... 6 more
Error: Not able to write product file: 'D:\S3\L2_WFR_subset\MANACAPURU\S3A_OL_2_WFR____20191104T135002_20191104T135302_20191124T134508_0179_051_124_3060_MAR_O_NT_002_subset.txt'
SEVERE: org.esa.snap.core.util.SystemUtils$SnapImagingListener: JAI error occurred: 'Problem occurs when computing a tile by the owner.' at com.sun.media.jai.util.SunTileScheduler@e17dd41
org.esa.snap.core.gpf.OperatorException: Not able to write product file: 'D:\S3\L2_WFR_subset\MANACAPURU\S3A_OL_2_WFR____20191104T135002_20191104T135302_20191124T134508_0179_051_124_3060_MAR_O_NT_002_subset.txt'
at org.esa.snap.core.gpf.common.WriteOp.computeTile(WriteOp.java:360)
at org.esa.snap.core.gpf.internal.OperatorImage.computeRect(OperatorImage.java:80)
at javax.media.jai.SourcelessOpImage.computeTile(SourcelessOpImage.java:137)
at com.sun.media.jai.util.SunTileScheduler.scheduleTile(SunTileScheduler.java:904)
at javax.media.jai.OpImage.getTile(OpImage.java:1129)
at com.sun.media.jai.util.RequestJob.compute(SunTileScheduler.java:247)
at com.sun.media.jai.util.WorkerThread.run(SunTileScheduler.java:468)
Caused by: java.lang.ArrayStoreException
at java.lang.System.arraycopy(Native Method)
at org.esa.snap.csv.dataio.writer.CsvProductWriter.writeBandRasterData(CsvProductWriter.java:142)
at org.esa.snap.core.gpf.common.WriteOp.writeTileRow(WriteOp.java:408)
at org.esa.snap.core.gpf.common.WriteOp.computeTile(WriteOp.java:326)
... 6 more
SEVERE: org.esa.snap.core.util.SystemUtils$SnapImagingListener: JAI error occurred: 'Problem occurs when computing a tile by the owner.' at com.sun.media.jai.util.SunTileScheduler@e17dd41
org.esa.snap.core.gpf.OperatorException: Not able to write product file: 'D:\S3\L2_WFR_subset\MANACAPURU\S3A_OL_2_WFR____20191104T135002_20191104T135302_20191124T134508_0179_051_124_3060_MAR_O_NT_002_subset.txt'
at org.esa.snap.core.gpf.common.WriteOp.computeTile(WriteOp.java:360)
at org.esa.snap.core.gpf.internal.OperatorImage.computeRect(OperatorImage.java:80)
at javax.media.jai.SourcelessOpImage.computeTile(SourcelessOpImage.java:137)
at com.sun.media.jai.util.SunTileScheduler.scheduleTile(SunTileScheduler.java:904)
at javax.media.jai.OpImage.getTile(OpImage.java:1129)
at com.sun.media.jai.util.RequestJob.compute(SunTileScheduler.java:247)
at com.sun.media.jai.util.WorkerThread.run(SunTileScheduler.java:468)
Caused by: java.lang.ArrayStoreException
at java.lang.System.arraycopy(Native Method)
at org.esa.snap.csv.dataio.writer.CsvProductWriter.writeBandRasterData(CsvProductWriter.java:142)
at org.esa.snap.core.gpf.common.WriteOp.writeTileRow(WriteOp.java:408)
at org.esa.snap.core.gpf.common.WriteOp.computeTile(WriteOp.java:326)
... 6 more
SEVERE: org.esa.snap.core.util.SystemUtils$SnapImagingListener: JAI error occurred: 'Problem occurs when computing a tile by the owner.' at com.sun.media.jai.util.SunTileScheduler@e17dd41
org.esa.snap.core.gpf.OperatorException: Not able to write product file: 'D:\S3\L2_WFR_subset\MANACAPURU\S3A_OL_2_WFR____20191104T135002_20191104T135302_20191124T134508_0179_051_124_3060_MAR_O_NT_002_subset.txt'
at org.esa.snap.core.gpf.common.WriteOp.computeTile(WriteOp.java:360)
at org.esa.snap.core.gpf.internal.OperatorImage.computeRect(OperatorImage.java:80)
at javax.media.jai.SourcelessOpImage.computeTile(SourcelessOpImage.java:137)
at com.sun.media.jai.util.SunTileScheduler.scheduleTile(SunTileScheduler.java:904)
at javax.media.jai.OpImage.getTile(OpImage.java:1129)
at com.sun.media.jai.util.RequestJob.compute(SunTileScheduler.java:247)
at com.sun.media.jai.util.WorkerThread.run(SunTileScheduler.java:468)
Caused by: java.lang.ArrayStoreException
at java.lang.System.arraycopy(Native Method)
at org.esa.snap.csv.dataio.writer.CsvProductWriter.writeBandRasterData(CsvProductWriter.java:142)
at org.esa.snap.core.gpf.common.WriteOp.writeTileRow(WriteOp.java:408)
at org.esa.snap.core.gpf.common.WriteOp.computeTile(WriteOp.java:326)
... 6 more
15/09/2020 17:02:00 - GPT processing complete, output file should be generated at:
D:\S3\L2_WFR_subset\MANACAPURU\S3A_OL_2_WFR____20191104T135002_20191104T135302_20191124T134508_0179_051_124_3060_MAR_O_NT_002_subset.txt
Finished in 15.44 second(s).
Unfortunately it crashed before extracting the pixel values and I end up with an empty text file (It happened both on Ubuntu18.04 and Windows 10 Pro).
I know this is no trivial matter so any help on the subject would be greatly appreciated.
These are the auxiliar files:
gpt_manacapuru.xml (1.1 KB) manacapuru.dbf (259 Bytes) points.prj (145 Bytes) manacapuru.shp (820 Bytes) manacapuru.shx (108 Bytes) manacapuru.kml (2.4 KB)
I can also try uploading the original S3 image somehow if you need it.
Thanks again,
David
You said it was working before, right?
Have you checked if that particular SEN3 is corrupted?
Another fast test would be if you could save your first subset in .dim and not .txt. But if it worked before I guess it could be more related to problems with specific images.
Hello and thanks for your suggestions @abruescas.
Yes, it was working before, I have also checked the xml graph file for any inconsistencies but I couldn’t find any )-:
I have also tested with several other images and they all look fine as you can see them opened in SNAP:
In the above image you may be able to see a greenish vector around the Manacapuru river in the Amazon. I’ve managed to extract the pixels before from several images using an automated version of the example graph XML that @marpet gave me, and then I managed to generate my time series below:
Now I have downloaded the all-time images available for Sentinel-3 over my region of interest to expand my series, but whenever I call the GPT tool either from inside my code or even manually, I always get the same error stack. I am honestly clueless…
It seems you have to import the shape file before subsetting. Send you the xml that worked for me. gpt_manacapuru_v2.xml (1.2 KB)
I saved the file in .dim.
I have been checking and the only way you can save your masked data in a txt file is using the Export/Other/Mask Pixels and I think is not available in gpt. If you just need a mean of the pixels in the ROI you could check the StatisticsOp.
