Support: 20/July/2007, Chris Hecker, WM2006/06

[mggr]

Contacted by Chris Hecker (also previously asked about azgcorr versions in #38). Full email below, but the key points are:

1. I want to use the ARSF/NERC LIDAR data that was recorded simultaneously with the hyperspectral data to improved geocorrection results. The data I want to correct are Eagle/Hawk data strips. The LIDAR data comes as point clouds in files like "Str_287.all". Can these point clouds be read directly into azgcorr (as what is referred to as an "external file") or do I need to use another software (e.g. SCOP++) to make regular gridded altitude data first?
2. For the LIDAR data, I presume that I wont need any geoid - ellipsoid correction parameter, since aircraft navigation data as well as airborne LIDAR are both based on DGPS heights above WGS84 ellipsoid. Is that correct?
3. I saw that there are a lot of special CASI options mentioned in the manual. Are there any for the SPECIM sensor in this version already and where would I find a list of these options?
4. Now here comes the tricky one:

   In the help file there is the option -E which allows the preprocessing of the LIDAR data cloud to new coordinates. Since all our basemaps are in European datum 1950, I wanted to pre-process the LIDAR cloud to this datum as well.
   Currently the data is in the standard format (UTM zone 30, WGS84) with xcoodrinate including the UTMzone as a prefix.
   The output datum is again UTM zone 30 with International ellipsoid and the following datum parameters dx=-87 dy=-98 dz=-121. no rotational parameters are known to me.
   
   So I tried the following AZGcorr statement:
   azgcorr -ELpt f1 x -d7 0 -87 -98 -121 0 0 0 0 -mUTMZ 30 -e Str_287.all
   
   Which to my knowledge should have given me a properly refernced data cloud for my UTM / ED50 system. I can see from the coordinates that are displayed during the processing that it is not going OK. The conversion of the x coordinate to the longitude value is already going wrong. I also noticed that the central meridian (in this case -0.0396 but changing with each time azgcorr is run ) doesnt make any sense. Am I doing something wrong or is this a bug? I am running azgcorr v. 4.8.3 on a 486pc-linux system.
   

Full email:

Hi Mike;
I started to work myself through the various manuals and LINUX version help file of azgcorr.
I have some questions regarding the processing. Maybe you can help me out with this. I am sure they (i.e. the first three) are pretty straight forward for you.

A) I want to use the ARSF/NERC LIDAR data that was recorded simultaneously with the hyperspectral data to improved geocorrection results. The data I want to correct are Eagle/Hawk data strips. The LIDAR data comes as point clouds in files like "Str_287.all". Can these point clouds be read directly into azgcorr (as what is referred to as an "external file") or do I need to use another software (e.g. SCOP++) to make regular gridded altitude data first?

B) For the LIDAR data, I presume that I wont need any geoid - ellipsoid correction parameter, since aircraft navigation data as well as airborne LIDAR are both based on DGPS heights above WGS84 ellipsoid. Is that correct?

C) I saw that there are a lot of special CASI options mentioned in the manual. Are there any for the SPECIM sensor in this version already and where would I find a list of these options?

D) Now here comes the tricky one:
In the help file there is the option -E which allows the preprocessing of the LIDAR data cloud to new coordinates. Since all our basemaps are in European datum 1950, I wanted to pre-process the LIDAR cloud to this datum as well.
Currently the data is in the standard format (UTM zone 30, WGS84) with xcoodrinate including the UTMzone as a prefix.
The output datum is again UTM zone 30 with International ellipsoid and the following datum parameters dx=-87 dy=-98 dz=-121. no rotational parameters are known to me.

So I tried the following AZGcorr statement:
azgcorr -ELpt f1 x -d7 0 -87 -98 -121 0 0 0 0 -mUTMZ 30 -e Str_287.all

Which to my knowledge should have given me a properly refernced data cloud for my UTM / ED50 system. I can see from the coordinates that are displayed during the processing that it is not going OK. The conversion of the x coordinate to the longitude value is already going wrong. I also noticed that the central meridian (in this case -0.0396 but changing with each time azgcorr is run ) doesnt make any sense. Am I doing something wrong or is this a bug? I am running azgcorr v. 4.8.3 on a 486pc-linux system.
Below I copy a screendump.

Thanks for the help,

Chris

His dump:

hecker@x7:~/in$ azgcorr -ELpt f1 x -d7 0 -87 -98 -121 0 0 0 0 -mUTMZ 30 -e Str_287.all
 
-----------------------------------------------------------------------
azgcorr  -- ver: 4.8.3-lin Jun 12 2007   (C) Azimuth Systems UK 1996, 2007
 

DEM pre processing...
 
Initialised transform for point cloud to WGS84 geogs...
 
Projection and datum shift details initialised....
Spheroid name: WGS84
  semi-major: 6378137.00000  semi-minor: 6356752.31425  e^2: 0.006694380 1/f: 29
8.257223563
Projection name: UTM
  origin  lat: 0.0000  long (central meridian) : -0.0396
  grid coords at origin  easting: 500000.00  northing: 0.00
  scale factor: 0.99960000  hemisphere: north
Datum shift name: none
 
Initialised transform for point cloud WGS84 geogs to local datum and projection.
..
 
Projection and datum shift details initialised....
Spheroid name: INTERN
  semi-major: 6378388.00000  semi-minor: 6356911.94613  e^2: 0.006722670 1/f: 29
7.000000000
Projection name: UTM
  origin  lat: 0.0000  long (central meridian) : -3.0000
  grid coords at origin  easting: 500000.00  northing: 0.00
  scale factor: 0.99960000  hemisphere: north
Datum shift name: SING
  old spheroid name: INTERN
  semi-major: 6378388.00000  semi-minor: 6356911.94613  e^2: 0.006722670 1/f: 29
7.000000000
  dx: -87.0000 dy: -98.0000 dz: -121.0000 metres  sc: 0.000000  ppm
  rx: 0.0000 ry: 0.0000 rz: 0.0000 secs
 

old x: 593792.77 y: 4095708.38  z: 131.270 lat: 37.0028547 lng: 1.0145276
new x: 857185.96 y: 4102748.87 z: -12.407
old x: 593792.77 y: 4095708.39  z: 131.310 lat: 37.0028548 lng: 1.0145276
new x: 857185.96 y: 4102748.88 z: -12.367
old x: 593792.79 y: 4095707.83  z: 131.270 lat: 37.0028497 lng: 1.0145278
new x: 857186.00 y: 4102748.32 z: -12.407
old x: 593792.79 y: 4095707.84  z: 131.290 lat: 37.0028498 lng: 1.0145278
new x: 857186.00 y: 4102748.33 z: -12.387
old x: 593792.82 y: 4095707.10  z: 131.300 lat: 37.0028432 lng: 1.0145280
new x: 857186.05 y: 4102747.59 z: -12.377
old x: 593792.82 y: 4095707.12  z: 131.340 lat: 37.0028433 lng: 1.0145280
new x: 857186.05 y: 4102747.61 z: -12.337
old x: 593792.91 y: 4095706.78  z: 131.310 lat: 37.0028403 lng: 1.0145290
new x: 857186.15 y: 4102747.27 z: -12.367
old x: 593792.90 y: 4095706.80  z: 131.370 lat: 37.0028404 lng: 1.0145289
new x: 857186.14 y: 4102747.29 z: -12.307
 
file: Str_287.all  converted points: 3910381
 
End of DEM run

[mggr]

Spoke to Gabriel to warn him about potential LIDAR queries, then cc'ed him on the reply. Also sent him a copy of azgcorr (he was interested in the LIDAR point cloud options) and asked him to let us know of any contact with Chris.

Replied to Chris on his questions 1-3:

1. needs gridded data, suggested GRASS, ARC or other GIS packages
2. yes @ WGS-84
3. not sure on the SPECIM options - suggested trying the CASI ones if needed (as the sensors are similar and it works this way in the internal processing). Offered to contact Bill if needed. Another potential contact would be Andrew.

On the 4th question (D), I said we'd try to reproduce it and get back to him. Sounds like a possible bug.

[mggr]

Reproduced the same thing here - I get a different central meridian (for the input datum/projection?) on each run.

Told Chris we'd contact Bill. In the meantime, suggested he used an external package for datum conversion and gridding.

[mggr]

Bug ticket #43

[mggr]

Hawk processed without a DEM

[mggr]

Hawk processed with a DEM

[mggr]

Followup query: Chris has had some weird problems when using a DEM.

Using `azgcorr -1 h153131b.hdf -3 h153131b_out_noDEM.hdf -BIout -dNO -mTM 3 0 0
-3 0.9996 0 500000 -in -p 5 5 -bl 100 -1`, he gets:

but when processing with a DEM that he's made (azgcorr -1 h153131b.hdf -3 h153131b_out_DEM2.hdf -BIout -dNO -mTM 3 0 0 -3 0.9996 0 500000 -in -p 5 5 -bl 100 -1 -e Str_287_all_fst_WGS84_02m.dat -es NO -ed 0 2675 468 593790 4095586 599138 4096520 2) he gets:

His query:

Some of the pixels seem OK but most of them are set to 0.0. In some
cases (if I use ED50 for DEM and output file) the northern half is also
entirely missing. Is that an effect of not using the -ef -it or -g flags?

I suggested he try these, but suspect that it's more likely to be a problem with the DEM. Asked him to re-run with -v while I got hold of the data locally.

[mggr]

Acquired Hawk data for this specific flight line, putting in ~/support/20070720-WM2006_06-_-ChrisHecker/hawk_blackout

[mggr]

Noticed that the DEM has height values below zero - this will probably cause problems with the -ez flag:

 -ez v      : DEM grid empty nodes value, default: 0.00
            :   NB: must be lower than all valid DEM node values

Also, while processing other flights, we noticed that a DEM that doesn't cover the full area will lead to black cut-outs.

Passed both recommendations to Chris.

[mggr]

mggr's processing of the Hawk with ASCII DEM

[mggr]

Chris came back with continuing problems and more diagnostics indicating the DEM wasn't read correctly. He also put the DEM on an FTP site (downloaded to ~arsf/support/20070720-WM2006_06-_-ChrisHecker/hawk_blackout/hecker_dem).

Looking it over, he has a binary DEM and it seems that azgcorr was trying to read it as an ASCII one (guessing). Converting the DEM to ASCII with ENVI (and editing a little) was enough to make things work with:

azgcorr -1 h153131b.hdf -3 h153131b_out_DEM2.hdf \
        -BIout -dNO -mTM 3 0 0 -3 0.9996 0 500000 \
        -in -p 5 5 -bl 100 -1 \
        -e Str_287_all_fst_WGS84_02m.txt -es NO \
        -ed 0 2675 468 593790 4095586 599138 4096520 2 -v 

(note the imagery stops at the edge of the DEM)

Looking deeper, I think he needs to use -eb and -edx rather than -e and -ed. Not got it working yet though, so emailed him about the ASCII variant.

[mggr]

Working version for binary DEM:

  azgcorr -1 h153131b.hdf -3 h153131b_out_DEM2-bin.hdf \
          -BIout -dNO -mTM 3 0 0 -3 0.9996 0 500000 -in \
          -p 5 5 -bl 100 -1 -v \
          -es NO -eb Str_287_all_fst_WGS84_02m.dat \
          -edx 0 2675 468 593790 4095586 599138 4096520 2 0 0 1 1.0 0.0 

[mggr]

Chris confirms he has the DEM stuff working, but now has run into problems with a datum shift:

Since all our datasets are in ED50, I was trying in a next step, to get
that working. However, azgcorr seems to have troubles with the -d7
option. I use the following statement (without DEM at the moment):

azgcorr -1 h153131b.hdf -3 h153131b_out_noDEM_ED50.hdf -BIout \
-d7 0 -87 -98 -121 0 0 0 0 -mUTMZ 30 \
-in -p 2 2 -bl 100 -1 \
-v > h153131b_out_noDEM_ED50.txt

This should give me UTM Zone 30 for ED50 mean. However, the -v output to
the txt file shows that it uses Datum Shift type "Dutch" which is made
for RD reference system (not UTM at all) with Bessel ellipsoid etc, so
something different.

Reproduced locally and don't understand it either ;) Tried various combinations and orders of options, but without much luck. Going to contact Bill on this one - ticket #44.

[mggr]

Contact from Chris Hecker again:

Hi Mike.
I was working on a paper last week, that's why you heard so little from me! ;-)
 
I am now back to processing the data. I have several questions for you. maybe you can figure out some of it.
 
a) had some black bars in data. They usually occur close to the periphery.  In fact, some of the pixels are shifted too far out by the topographic correction.  The pixels beyond the black line are not in the correct place. They should be adjacent to the rest of the strip (see attached jpg).  I can bridge the gap by setting -g to higher values, but that is not what I want. I would like to have them in the correct place! ;-) I tried to limit the output to a fixed size grid, process only the input lines that fall into the DEM (DEM strip is shorter than hawk strip) and increased the space around the dem (filled with 0.0) in order to change it.  Nothing helped except that increasing the DEM spacing changed the location of the black strip.
 
Just a few minutes ago i found out that I had set the nodata value to -999 while the DEM still contained 0.0 values for NoData. When I mask them and set to 999.0 the black line is gone. The pixels beyond the black line have been totally removed though. However, this means that AZGCORR, while calculating the position of a pixel that is still WITHIN the extent of the DEM, is thrown off by an altitude value (0.0) that is beyond the intersection of view vector and topography.
There are two reasons I can think of but both dont make me happy:
It could be that it's using the second intersection of the view vector with the topography or
that there is misalignment between DEM and hawk data.
The second option is quite possibly the cause, since the pixels that used to be south of the black line, are now totally gone (possibly since azgcorr thinks that they are in the part of the DEM that contains NoData values.
Puuuhhh .....
 
Here is one of the CLs I used:
azgcorr -BIout -dNO -mTM 3 0 0 -3 0.9996 0 500000 \
-in -p 2.5 2.5 \
-bl 100 -1 \
-l 950 2850 \
-1 h153131b.hdf -3 h153131b_WGS84_bnd100.hdf \
-eb Str_287_all_fst_WGS84_halfm.dat \
-edx 0 10700 1870 593790 4095585.5 599139.5 4096520 0.5 0 1 1 1 -999.0 \
-es NO \
-v > h153131b_WGS84_bnd100.txt
 
I put the files I was using onto the ftp server. Maybe that helps in re-producing.

ftp.itc.nl/pub/hecker

Username: hecker

Password: BmpBze
 
 
b) I have a question about the -u functions of adjusting the NAV data. Is this still necessary? I have noticed small absolute geolocation discrepancies in 1 flightline but I am not sure if it is constant through all scenes.  In would expect that these small corrections would be identical for all flightlines recorded right after each other (or at least that they slowly change in case they are related to the drift of the IMU).  Is that assumption correct?  If yes, it would of course save a lot of time as I dont have to figure out the correction values for each strip.
 
c) I understand that the LIDAR point data were processed by somebody else. I am missing the LIDAR data for one of our flightlines. I am not sure, whether they were not recorded or whether they were missed in the processing. Is there a list somewhere that shows which instruments were on for each flightline?  The flightline (in hawk "codes") is flightline h153121b.hdf of the Mojacar study area. Its a pretty short line.
 
d) I had quite a hard time cleaning the LIDAR ascii data up since the data contained empty lines, some binary codes (looks like some ?C+ codes slipped into the ASCII file) etc. I am not sure if that is maybe a corruption that came during ftp, but maybe they should know about this.
 
e) here is another bug report:
the extended DEM header option - edx doesnt work with re=1.
even though the number of columns and lines are correct for the coordinates, the message :

DEM file: Str_287_all_fst_WGS84_02m_undef.dat
    limits SW: 593789.0 4095585.0  NE: 599139.0 4096521.0
    order (0): S-N cols: 2675 rows: 468 grid inc: 2.0 fmt: 1
** info: one or more grid limits are NOT multiples of the grid increment

is displayed. seems that the re=1 flag is not read properly.
I also think that the -edx or=1 flag is not read properly but I havent tested that in detail.
maybe you can forward this to the programmer. 

Followed by:

just noticed that I made a big mistake in the CL. I think I had the -edx or = 0 while in ENVI it has to be or=1. So it read my DEM upside down!  That's probably why the pixels ended up being outside the area covered by the DEM although they werent really close to the edge.
 
so entire a) point of mine as well as the second half of e) (that -edx or = 1 may not be read properly)  can be deleted! ;-)

Points that remain:

• geolocation corrections
• lidar missing
• corruption in the lidar imagery
• DEM limits warning

Sent responses on all these as follows:

Hi Chris,

Good to hear you're back at the data ;)  Bill, the software author, is back from Iceland too so we're just following up on your earlier queries - I'll let you know how those go.

Chris Hecker wrote:
> b) I have a question about the -u functions of adjusting the NAV data. Is this still necessary? I have noticed small absolute geolocation discrepancies in 1 flightline but I am not sure if it is constant through all scenes.  In would expect that these small corrections would be identical for all flightlines recorded right after each other (or at least that they slowly change in case they are related to the drift of the IMU).  Is that assumption correct?  If yes, it would of course save a lot of time as I dont have to figure out the correction values for each strip.

With regard to the geolocation, the normal approach at ARSF is to geolocate the strips as accurately as we can with what we have.  If the flights are in the UK, we normally use the Ordnance Survey vector overlays and a DEM to ensure we have good geolocation.  In flights abroad, we often don't have access to vector overlays or DEMs (at least until the LIDAR is available anyway), so have to rely on previous calibrations and visual checks - flights like this will probably be less accurately geolocated, but we can't quantify.

Typically, all the flightlines will have the same residual error, so you should be fine doing one correction and applying to all.

The only exceptions to this would be if there was something like a system crash necessitating reboots of the sensors or similar problems. I believe that some of the Specim sensors would have problems of synchronisation following dropped scanlines, but I'm not sure if these would have affected your data (I can make enquiries if you wish?).

You'd normally be able to notice any problems like the above as you'd find that straight lines on the ground will suddenly go wobbly-looking.  Quickly looking over your data should be good enough to determine if you have this problem or not.  I'd expect you're probably fine to use a single correction.

> c) I understand that the LIDAR point data were processed by somebody else. I am missing the LIDAR data for one of our flightlines. I am not sure, whether they were not recorded or whether they were missed in the processing. Is there a list somewhere that shows which instruments were on for each flightline?  The flightline (in hawk "codes") is flightline h153121b.hdf of the Mojacar study area. Its a pretty short line.

LIDAR is processed by the Unit for Landscape Modelling at the Uni. of Cambridge (http://www.uflm.cam.ac.uk/).  They have a data processing status page (http://www.uflm.cam.ac.uk/dataproc-completed.html) that lists your project as complete.

It also has a thumbnail of the final DEM at http://www.uflm.cam.ac.uk/thumbnails%5C06153b_Mojacar.jpg - perhaps you can spot if they have something you don't?

I've also attached the flight logsheet which lists what was flown.  If something was flown and the LIDAR appears to be missing, the best thing would be to contact ULM to see if they know what happened.  The contact person there is Gabriel Amable (gsa1000@cam.ac.uk).

> d) I had quite a hard time cleaning the LIDAR ascii data up since the data contained empty lines, some binary codes (looks like some ?C+ codes slipped into the ASCII file) etc. I am not sure if that is maybe a corruption that came during ftp, but maybe they should know about this.

Hmm, that's definitely interesting - again, I'd suggest contacting Gabriel, but would appreciate it if you'd cc us in on it (or at least the outcome) so we know what to look for.

> e) here is another bug report:
> the extended DEM header option - edx doesnt work with re=1.
> even though the number of columns and lines are correct for the coordinates, the message :
>
> DEM file: Str_287_all_fst_WGS84_02m_undef.dat
>     limits SW: 593789.0 4095585.0  NE: 599139.0 4096521.0
>     order (0): S-N cols: 2675 rows: 468 grid inc: 2.0 fmt: 1
> ** info: one or more grid limits are NOT multiples of the grid increment

I think the warning here relates to the grid limits not being an integer multiple of the grid increment (ie. they should be divisible by 2).  One would think that the fact the range is divisible would be enough though..

Hope this helps - keep us updated on how you do :)

Also, I'm likely to take a couple of days holiday at the end of this week, so if you email me during that time, I might take a little while to reply.  It'd probably be best if you cc arsf-processing@pml.ac.uk on all emails - that's an email alias that goes to about 6 people.  That way you'll be sure that someone at least is reading any urgent queries :)

Cheers,

Mike.

[mggr]

Sent (Wed 15th Aug) the latest azgcorr to Chris, which may fix the problems he had. He's busy atm but will respond when he gets time.

[mggr]

He's working on it a bit more and phoned with a couple of questions that are Bill-level:

• First off, he was running with a 50cm DEM (about 20kx40k pixels!) and running into memory problems, so I advised he restrict that to 2m accuracy at most (which is the pixel size he was aiming for). He also said he was getting an "error 6" when splitting his DEM into two files, but I'm not sure what that is ;)
• The second question was what happens if his DEM is offset slightly from the grid - it gives a warning about this, but does it interpolate?

Sending Bill an email asking if he's ok being contacted directly.

[mggr]

Some communication with Bill, but gone quiet now (possibly due to fieldwork).

Contact Chris in a little while to see if he's happy.

[mggr]

New contact in #108. Closing this ticket as he seems happy with the data.