dongjief / ledaps

Automatically exported from code.google.com/p/ledaps
2 stars 1 forks source link

Discrepancy in Ledaps results #7

Closed GoogleCodeExporter closed 9 years ago

GoogleCodeExporter commented 9 years ago
I am having an issue wrt the results of Ledaps generated by my build of 1.2.1 
and the products ordered through  https://espa.cr.usgs.gov.
The values differ (for some images) in bands 1, 2,3 and 4 (only bands I have 
looked at) by between 5 and 100 which seems to be beyond what I should be 
expecting due to differences in architecture between machines. Should I be 
concerned with this? In most case the difference is less than 1%.

I was wondering if this might be due to differences in the ancillary files 
being used? I downloaded my ancillary files from 
http://landsat.usgs.gov/espa/files/. 

I started to investigate and attempt to build the files within the directory 
ledapsAncSrc. Is this useful or necessary? If no please ignore the following 
paragraph. :)

I am having issues as there appears to be no documentation to supplement this 
package. I tried to run the makefile, but failed with the message 
"convert_ozone.c:5:25: error: hdf4_netcdf.h: No such file or directory". A 
search of my system and the dependent libraries (installed for ledaps) for 
"hdf4_netcdf.h" was unsuccessful. I did however find a "netcdf.h" file and 
changed the include statements to netcdf.h instead in both the "convert_ozone.c 
and nep_reanalysis_...c as well as including the szlib in the library and 
include declarations in the makefile. And volia it compiled and installed, and 
the executables were copied into the ledapsSrc/bin file as expected.

I reran ledaps with the new executables - from ledapsAncSrc in place - but the 
results still did not jive with the NASA generated images. Even stranger still 
is the fact that the SR result of ledaps for LT50330242010197EDC00 matched the 
expected results within 1 BUT for LT50360232010218EDC00 the results are as 
stated above. Confounding as both images are from the same sensor and year. 

Can you give me any suggestions or insight as to why this discrepancy is 
occurring and whether I should be concerned?

As an aside. I should note that when I built ledaps I needed to edit several of 
the C files (ar.c, lndsr.c, myhdf.c, output.c and run_sixs.c) and expunge all 
the old school // commented code ( as opposed to the /* comment */ style) since 
my compiler would not accept these // commented lines.

I look forward to your response. Thanks for your assistance.

Original issue reported on code.google.com by kayak2bl...@gmail.com on 26 Jun 2013 at 6:28

GoogleCodeExporter commented 9 years ago
The ledapsAncSrc includes two scripts (updatencep.py and updatetoms.py) which 
allow the user to process and update the auxillary data files needed as input 
to the LEDAPS application.  The initial tarball of auxillary data files we have 
hosted on ESPA are not up-to-date with the most recent auxillary products.  
Thus, if you desire to process Landsat products newer than the last date in 
2012 from the ledaps_anc_1980-2012.tar.gz file, then you'll need to compile and 
run the scripts in the ledapsAncSrc to update your auxillary data files for 
ozone, pressure, water, etc.  If you only desire to process products between 
1984 and 2012, then you should be ok without running the auxillary update 
scripts.  However, there is a small chance that the underlying TOMS and/or NCEP 
files hosted by NOAA and NASA may be updated.  In that case, running the 
auxillary scripts will also update those files.

To process all the auxillary data back to the start of the data availability, 
the user would want to run the updatencep and updatetoms scripts using the 
--quarterly flag.  This is what the ESPA team runs on a quarterly basis.  If 
the user only desires to process/update auxillary data for the most recent 
year, then the user would want to use the --today flag.  And, if there are only 
particular years which the user desires, then the --start_year and --end_year 
flags should be used for the update scripts.  To update the auxillary data for 
LEDAPS, the user will want to run both updatencep.py and updatetoms.py.  The 
order of which script is run first does not matter.

The hdf4_netcdf.h file should be in the include subdirectory of your HDF4 
installation.

Regarding the differences in some of your outputs vs. those ordered from ESPA, 
I agree that some of the larger differences are likely not due simply to 
differences in architecture or operating systems.  I would recommend making 
sure to compare the latest version of ESPA results with the latest version of 
the source software.  There is actually a global attribute in the output 
lndcal, lndth, and lndsr HDF files which depicts the version of LEDAPS used to 
process the data.  (Ex.  LEDAPSVersion = 1.2.1).  This version number is 
identified in the lndpm.c file under src/lndpm.

With respect to the old-style comments, I'll have to look at getting all of 
those changed if they are creating problems for some users, such as yourself.  
Thanks for the heads-up.

Original comment by schmidtg...@gmail.com on 28 Jun 2013 at 1:41

GoogleCodeExporter commented 9 years ago
The latest version of the LEDAPS software now has the C++ style comments (//) 
replaced with C style comments (/* */).  This will be released with version 
1.2.3 in late July.

Original comment by schmidtg...@gmail.com on 28 Jun 2013 at 2:19

GoogleCodeExporter commented 9 years ago
[deleted comment]
GoogleCodeExporter commented 9 years ago
Hi Gail. Thanks for your response. WRT the discrepancy in the results between 
the USGS product and the product I generated with my build both were produced 
with version 1.2.1.  As I mentioned previously the curious thing is that SR 
result of ledaps for LT50330242010197EDC00 matched the expected results within 
1 digit BUT for LT50360232010218EDC00 the results were between 5 and 100 digits 
(around 1%) of the expected (USGS generated) values. Any ideas as to why this 
might be occurring? Thanks

Original comment by kayak2bl...@gmail.com on 2 Jul 2013 at 3:29

GoogleCodeExporter commented 9 years ago
Let me run some tests on my end with my built version of the software vs. the 
ESPA-generated products.

Original comment by schmidtg...@gmail.com on 2 Jul 2013 at 4:48

GoogleCodeExporter commented 9 years ago
Thank you Gail. I would recommend starting with the two images that I listed 
above as it seems to be a transient issue depending on the images used.

Original comment by kayak2bl...@gmail.com on 2 Jul 2013 at 6:41

GoogleCodeExporter commented 9 years ago
Hi Gail.
I have some other information that might be useful wrt to this issue.
I was doing some comparative analysis on cloud free overlapping forested areas 
between various images that were generated with my build and the USGS generated 
product. I found a common area between 3 LT5 images from different path rows 
that were cloud free comprised of 1715 pixels. I then ran up some basic stats 
for the extracted values for all the bands of the .hdf product.

I have a attached a small extract of my analysis that focuses solely on the AOT 
values. Two works sheets show the extracted values for all 1715 pixels and are 
named to reflect their source (Mine - Canadian Forest Service (CFS))ad the USGS 
product.
There is a Summary worksheet that compares the summary stats between these 
sources. 

Some of the he differences between the images generated by each processes is 
probably due to re-sampling and geometry, but what really jumped out was the 
difference between my AOT generated values and the USGS AOT values for the same 
images. As I understand it, this layer is generated from the LTI1 Band 1 
product so somehow the algorithm or some values within it are different between 
my build of 1.2.1 and the USGS 1.2.1 build.

I hope this helps and I also hope that you had a great time on the 4th of July.
Best regards. If you want to call and discuss any of this I am here from 8 - 4 
PDT today and next week from 5:30-3:00 PST. Thank again for your attention to 
this.

Geordie Hobart, BEng
Research Assistant
Natural Resources Canada / Ressources naturelles Canada
Canadian Forest Service / Service canadien des forêts
Pacific Forestry Centre / Centre de foresterie du Pacifique
506 W. Burnside Road / 506 rue Burnside Ouest
Victoria, BC  V8Z 1M5 / Victoria, C-B V8Z 1M5

Tel: (250) 298-2403    Fax: (250) 363-0775
mailto:ghobart@nrcan.gc.ca

Original comment by kayak2bl...@gmail.com on 5 Jul 2013 at 4:33

Attachments:

GoogleCodeExporter commented 9 years ago
Hi Geordie.  I processed your two scenes and was able to duplicate the results 
you are seeing, using the original auxiliary datasets that we provided on the 
LEDAPS website.  However, when I process both these scenes using the updated 
auxiliary data (via the LedapsAncSrc scripts), the differences for both scenes 
between the ESPA version are very, very minimal (a handful of pixels with 
differences of 1 for the entire lndsr product).  The lndth and lndcal products 
are exactly the same for ESPA and my local processing, using both the updated 
and original auxiliary data.  So, I believe you have correctly installed and 
built your software.  And the differences you are seeing between your 
processing and ESPA results are due to the differences in the auxiliary data.

We completely reprocess the auxiliary data every quarter in order to use the 
latest data available from NCEP and NASA.  It is our belief that if the file 
has been updated by those agencies, then we should use the latest and greatest 
data.  We also process the current days of data every night to make sure we 
have as up-to-date auxiliary data as possible for processing the most current 
Landsat products.

Having said all of this, I've done some further digging to see what kind of 
differences existed (in the area of each scene) within the ozone data.  
Unfortunately what I'm seeing is similar minor differences in the ozone data 
(values are off by one or equal) for BOTH scenes.  At this time, I believe 
there is an issue with how the ozone data is read and interpolated for 
processing through lndsr.  Stay tuned for more information and most likely a 
new release.

Original comment by schmidtg...@gmail.com on 9 Jul 2013 at 5:39

GoogleCodeExporter commented 9 years ago
the ozone data has changed over the years and from instrument to instrument.  
The pre-OMI products had a lower longitudinal resolution as compared to the OMI 
products.  The convert_ozone scripts are working correctly in building the HDF 
files that are then read into the LEDAPS code, both for the old and new 
resolutions.  Therefore the HDF files representing the auxiliary data is 
accurate in the representation of this higher resolution.

However, I have confirmed that the ozone metadata is hard-coded in the lndsr 
application (including the min/max lat/long values and the delta lat/long 
values).  The hard-coded values are for the older ozone data products, thus 
using the 1.25 degree delta from longitude to longitude.  Therefore, the new 
ozone data is not being read correctly since the min/max and and particularly 
the delta long values are different from the old values.

Even more interesting is that the ozone data is stored from S latitude to N 
latitude in the original text files from GSFC, and this is being flip-flopped 
in the auxiliary HDF files.  However, this is not being handled correctly when 
writing the latitude values to the HDF file.  The latitude values are still 
written as -89.5 to 89.5 from top of the scene to bottom of the scene.  For a 
while, I thought the ozone data was being read upside down, but in all 
actuality it's being read correctly and the latitude values are incorrectly 
represented.  I will get this fixed.

Here's an example of the various file headers for each of the instruments 
showing the min/max and delta values for the lat and longs.

 Day: 287 Oct 14, 1994 METEOR-3/TOMS ARC OZONE    GEN:05.074 V8 ALECT: 12:41 PM 
 Longitudes:  288 bins centered on 179.375 W to 179.375 E  (1.25 degree steps)  
 Latitudes :  180 bins centered on  89.5   S to  89.5   N  (1.00 degree steps)  

 Day: 365 Dec 31, 1978 NIMBUS-7/TOMS NRT OZONE    GEN:04.119 V8 ALECT: 11:50 AM 
 Longitudes:  288 bins centered on 179.375 W to 179.375 E  (1.25 degree steps)  
 Latitudes :  180 bins centered on  89.5   S to  89.5   N  (1.00 degree steps)  

 Day: 126 May  6, 1993 NIMBUS-7/TOMS NRT OZONE    GEN:04.119 V8 ALECT: 10:19 AM 
 Longitudes:  288 bins centered on 179.375 W to 179.375 E  (1.25 degree steps)  
 Latitudes :  180 bins centered on  89.5   S to  89.5   N  (1.00 degree steps)  

 Day: 316 Nov 11, 2004    OMI TO3    STD OZONE    GEN:12:096 Asc LECT: 01:47 pm 
 Longitudes:  360 bins centered on 179.5  W  to 179.5  E   (1.00 degree steps)  
 Latitudes :  180 bins centered on  89.5  S  to  89.5  N   (1.00 degree steps)  

 Day: 198 Jul 17, 2010    OMI TO3    STD OZONE    GEN:12:097 Asc LECT: 01:46 pm 
 Longitudes:  360 bins centered on 179.5  W  to 179.5  E   (1.00 degree steps)  
 Latitudes :  180 bins centered on  89.5  S  to  89.5  N   (1.00 degree steps)  

 Day: 218 Aug  6, 2010    OMI TO3    STD OZONE    GEN:12:097 Asc LECT: 01:48 pm 
 Longitudes:  360 bins centered on 179.5  W  to 179.5  E   (1.00 degree steps)  
 Latitudes :  180 bins centered on  89.5  S  to  89.5  N   (1.00 degree steps)  

So, I need to do two things.  First, I have to fix/update the lndsr code to 
handle the old ozone and new ozone resolutions.  Second, I have to update the 
convert_ozone code to correctly represent the fact that the data has been 
flipped and therefore the northern-most latitudes should be on the top.  And, I 
have to make sure I don't screw up the interpolation of the water vapor, 
surface pressure, and air temp, since they are processed through the same 
interpolation function. ;-)

One change is cosmetic (correctly representing the lat values in the auxiliary 
HDF files), and the other will allow the OMI ozone data to be processed 
correctly.  Unfortunately any data products from 2004 onward (when OMI was 
being used) processed through LEDAPS will have been incorrectly reading the 
ozone data (using incorrect min/max longitudes and incorrect longitudinal 
deltas).

Original comment by schmidtg...@gmail.com on 9 Jul 2013 at 7:20

GoogleCodeExporter commented 9 years ago
Thanks for the update Gail. Very impressive sleuthing on your part to track all 
this down. Good luck with the code updates! All the best.

Original comment by kayak2bl...@gmail.com on 10 Jul 2013 at 9:37

GoogleCodeExporter commented 9 years ago
A new version of the LEDAPS source code has been checked in.  It handles the 
discrepancies in the resolutions of the ozone data for OMI vs. pre-OMI 
platforms.

The auxiliary scripts have also been updated to correct the fact that the 
latitude values are flipped upside down and don't correctly represent the way 
the ozone data is written.  Again, this is more of a cosmetic change than 
anything, as LEDAPS will handle the flip-flopped values as well as the correct 
latitude values.

Please note that when processing using the new version there will be 
discrepencies between the new outputs and any previous outputs for products 
processed using the OMI ozone data (approximately 2004 onward).

The code has been released to the trunk as well as released to version 1.3.0.  
This version will be updated in ESPA at the end of July.

Original comment by schmidtg...@gmail.com on 12 Jul 2013 at 2:43

GoogleCodeExporter commented 9 years ago

Original comment by schmidtg...@gmail.com on 12 Jul 2013 at 2:43

GoogleCodeExporter commented 9 years ago
Gail. Thanks very much for all your efforts ... much appreciated. All the best.

Original comment by kayak2bl...@gmail.com on 12 Jul 2013 at 8:26

GoogleCodeExporter commented 9 years ago
Going to close this issue at this time.

Original comment by schmidtg...@gmail.com on 23 Jul 2013 at 6:30