/weathersat

HRPT and MODIS (image) data analysis s/w

Primary LanguageFortranGNU General Public License v3.0GPL-3.0

This directory contains the source files for HRPT and (AQUA) MODIS data 
processing, using a map-/geo-projection approach.

Also contained are the programs to process the .bin files of the X-band 
instruments (MERSI) on FY3B, FY3D and FY3E. These programs attempt to 
do:
- Data reading
- Fixed pattern removal
- Mirror side correction
- (pseudo) truecolour image generation

(There also is a beta version of some tracking s/w using what I call 
look-ahead, it specifically targets a SPID RAS/HR rotor with an MD-02 or 
MD-01 controller)

If you don't read this README with attention, as well as the 

./hrpt.exe --help 

output, you will (!!) fail to successfully run the s/w. Especially 
the environment variables described below are crucial !!!!

For the MERSI data processing there are specific instructions at the end 
of this file.

First and most important - it is not at all difficult to make the 
software crash - I tried to make it robust, but it has many options 
and by far not all are mutually compatible. Also, there are plenty of
bugs left, but it CAN work.

This software needs a UNIX environment. Either vintage booted Linux 
(possibly as a VM), or WSL under windows. Cygwin under windows is 
possible but does not support the large memory model and would thus 
not enable processing (AQUA) Modis imagery.
For WSL installation instructions : see below.

Inside the UNIX environment you need gcc, gfortran, ImageMagick and gmt.
This can be installed using (Ubuntu):

sudo apt-get install gcc gfortran imagemagick gmt

If your windows machine has < 16 GByte of memory you might not be
able to use the 'large' memory model (see below) - if that is the case 
I do not know if the s/w will be able to process 250 m MODIS data.
For all HRPT files a MEDIUM memory model is OK.

Currently the following are supported:

              From Oleg's tool       From XHRPT    From J-L deframed Modis
NOAA                                   .raw16
FY               .C10                  .bin
Metop            .hpt                  .bin
MN2              .dat                  .bin
Aqua MODIS                                              .bin
Terra MODIS                                             .bin

I know that for FY using .bin files has a bug that makes 1 in 15 files or so fail.
If so, use Oleg's s/w to convert to .C10 and have my s/w use the .C10 file.

Use ./hrpt.exe --help for command line explanation - this is essential to do !

To start:

Extract this GITHUB repository in full to a directory named HRPT.

Then get the Google drive file:

https://drive.google.com/file/d/1DV9jHR8bmoTfD5vDIGD6iiShFGfNFrrU/view?usp=sharing

and untar-gzip in the directory containing buildhrpt.sh. This will create a 
sub-directory resource/ with support files.

Then build the s/w, which is done as follows:

./buildhrpt.sh

Note this requires gfortran and running the s/w requires imagemagick (see WSL comments below).
This builds hrpt.exe and the support s/w polarstereo.exe, readbin.exe and readbin_modis.exe

The software requires 5 environment variables to run. In my WSL environment they are set 
as follows:

export TLE=/mnt/y/TLE/
export HRPTOUT=/mnt/f/Output/
export MMODEL=LARGE
export LAT=52.338989
export LONG=4.709772

set the MMODEL to MEDIUM if AQUA Modis processing is not required.
adapt LONG and LAT to your position - the setting above puts you smack on top of a Schiphol runway
for LONG - positive values are Eastern longitude

As the software is TLE based it needs a TLE file. I download these daily from 
the web site with a crontab script (you need to modify the local directories in here !) :

#!/bin/bash
mytlefile=/share/Public/TLE/`(date +%Y%m%d)`-weather.tle
wget http://www.celestrak.com/NORAD/elements/weather.txt
mv weather.txt $mytlefile
mytlefile=/share/Public/TLE/`(date +%Y%m%d)`-resource.tle
wget http://www.celestrak.com/NORAD/elements/resource.txt
mv resource.txt $mytlefile
mytlefile=/share/Public/TLE/`(date +%Y%m%d)`-geo.tle
wget http://www.celestrak.com/NORAD/elements/geo.txt
mv geo.txt $mytlefile

this creates the files:

20210306-weather.tle
20210306-resource.tle
20210306-geo.tle

my s/w uses the data from the HRPT filename to find the equivalent date TLE 
file (and also tries 1 day earlier and later) in that directory.

If you have set the TLE environment variable, this script could look like:

#!/bin/bash
mytlefile=$TLE`(date +%Y%m%d)`-weather.tle
wget --secure-protocol=auto https://www.celestrak.com/NORAD/elements/weather.txt
mv weather.txt $mytlefile
mytlefile=$TLE`(date +%Y%m%d)`-resource.tle
wget --secure-protocol=auto https://www.celestrak.com/NORAD/elements/resource.txt
mv resource.txt $mytlefile


This directory is defined by the TLE environment variable defined as above.
The output from the programs goes to HRPTOUT, see above for my default.


The hrpt.exe program can generate 

- Per channel B/W images
- A merged 3 channel image and a simple geometry corrected version
- A polar stereographic projected 3 channel colour version
- A few images indicating the swath of the observation and the subsatellite tracked
  (with an option to also show a defined solar zenith angle 'sunset' line)
- An intensity histogram plot per channel
- A cumulative solar zenith angle distribution, used when defining a merge with an 
  IR channel
- Merged visual and IR channel images

Usually the program runs as:

./hrpt.exe filename 221 gapcor project=221 

you can add a border and skip some not essential output:

./hrpt.exe filename 221 gapcor project=221 border fast

or also a gamma:

./hrpt.exe filename 221 gapcor project=221 border fast gamma=1.2

For the FY and NOAA-19 satellites the borders align well with the images. For 
some reason for NOAA-18 and Metop this is less good. For this there is a 
theta= option. For Southern passes this is negative, for northern passes 
positive. This value represents a rotation of the scan angle (in degrees) 
around the nadir vector to the subsatellite point.

I have tried to set defaults for all parameters which are reasonable - if only 
as a starting point to experiment. They are indicated when you run with the --help.

The szalim, szamerge and szafact define the IR channel merging. Example:

./hrpt.exe filename 221 gapcor project=221 gamma=0.95 szalim=83.0 szamerge=4 szafact=1.15

(you can off course also switch on the border and theta= command)
This would merge the channel 4 IR image multiplied by a factor 1.2 (on top of 
an auto-normalisation) for all pixels where the solar zenith angle is above 85 
degrees. There actually is a 'soft' gradual merge from a few degrees below 
to a few degrees above the user defined value.

The merge option allows to use the first argument as a text file listing the files to be 
merged. Like:

./hrpt.exe ./examples/mergefiles.txt 221 project=221 gamma=0.95 histcor merge

where mergefiles.txt contains (from an UBUNTU Oracle virtualbox example) :

../Win7/HRPT-Data/Decoded/FY/FY3C_2019-01-19_0945.C10 theta=0.01 linedelay=-6
../Win7/HRPT-Data/Decoded/Metop/2019-01-19_0959_M02.hpt theta=-2.7 linedelay=3

the theta values are what you derive from the hrpt program as outlined 
above. IMPORTANT: If you define theta for one file it HAS to appear for all files.
If you want to use theta as 0 degrees, specify e.g. theta=0.01 but NOT theta=0.0 
as this has a special meaning inside the s/w (an experimental latitude based auto 
correction).
The text file defining which passes to merge can contain five lines at most !

The linedelay option is a subtlety. If on the left side of the images the borders 
are below the country in the image and on the right hand side the other way around, 
you need theta to correct. If there is an apparent left and right symmetrical 
displacement of the borders with the actual countries, you can use linedelay 
to change the timestamp of all scanlines by a multiple of 167 msec (different for 
AQUA MODIS). The FY satellites sometimes suffer from this. A small value should do.

An example of a MODIS processing command line, using files from Jean-Luc:

./hrpt.exe ./examples/panmodis.txt 221 gapcor project histcor modispan=correct binary merge reflong=70.0 border=highres

with panmodis.txt containing :

/mnt/f/AQ_MODIS-2019-09-29-1700.bin
/mnt/f/AQ_MODIS-2019-09-29-1700.bin
/mnt/f/AQ_MODIS-2019-09-29-1700.bin

This somewhat silly construct is because it uses the batch/merge concept that already 
was built in to the S/W.

I will try to answer questions you have. Also, most likely I will need a copy of your 
file to try and reproduce the problem.

Finally, the s/w is reasonably fast (a few minutes per HRPT file), but an AQUA MODIS 
true colour with the modispan option easily takes 45 minutes and if you add sharpen, 
add another hour (but its worth it !!)

In WSL2 the s/w runs significantly faster than in WSL1.

TLETRACK Software:

The tletrack.exe is a dedicated program that allows using the ROT2PROG protocol through 
a comport (default COM4 if you need a different one in the s/w search for ttyS4 , replace 
the 4 with the COM port number you need and rebuild). It is specific for the MD-02 
controller of a SPID RAS/HR or a SPID BIG-RAS/HR rotor.

./tletrack --help 

will show the help text on how to use the software. It also allows to use JPL Horizons 
based tracking (with the simplification that it uses 1 minute tracking intervals, 
but for a dish say up to 2-3 m this is OK). If you use the JPL Horizons based tracking 
the lookahead feature currently is not useable.



DEPROJECT Software:

The deproject software allows to correct a raw RGB image (png or jpg) for Earth 
curvature AND - for AQUA and TERRA MODIS and FY3B MERSI-1, FY3D MERSI-2 and 
FY3E MERSI-LL to perform the Bow-Tie (approximate but pretty good) correction.

./deproject --help

will show the help text on how to use the s/w. An example of correcting an AQUA 
MODIS image is:

./deproject.exe /mnt/y/Test/MODIS-RGB-121.png alt=708.4 fov=110.0 bowtie=mersi1 bowgamma=1.6 theta1=0.0 theta2=1.05

yes, I know the bowtie=mersi1 here looks silly, but it enables the bowtie algorithm 
and the 3 parameters that follow set it to the AQUA/TERRA MODIS values.
For FY3B and 3D you can remove these 3 parameters and use the default.

NOTE: The s/w is tested and validated on both FY and MODIS high resolution channel 
RGB (!) images whose # of image lines is (and HAS TO BE) an integer multiple of 
40 lines. It should in principle also work for medium res and low res images IF you 
add the relevant parameter (see the help text).
If you just wish to use the earth curvature correction the 40 lines multiple requirement 
does not apply.


FY3xREAD Software:

The 3 programs fy3bread.exe, fy3dread.exe and fy3eread.exe expect .bin files as produced 
by the s/w of @aang254 (Twitter) or @LucMillette (Twitter).
These programs internally call the deproject software to apply earth curvature and 
Bow-Tie corrections.

The usage for the miscalaneous missions would be:

./fy3bread.exe /mnt/c/FY3B/mersi.bin correct

or:

./fy3dread.exe /mnt/c/FY3D/mersi.bin correct

or:

./fy3eread.exe /mnt/c/FY3E/mersi.bin correct invert

The invert switch is needed for FY3E to apply the fixed line pattern corrections in a 
(reasonably) successfull way. The output of these programs can be found in the directory 
where the input file resides.



TROUBLESHOOTING : 

In case of imagemagick complaints about limits being exceeded:

cat /etc/ImageMagick-6/policy.xml

and check the constraints. I now use in WSL :

  <policy domain="resource" name="memory" value="1GiB"/>
  <policy domain="resource" name="map" value="2GiB"/>
  <policy domain="resource" name="width" value="32KP"/>
  <policy domain="resource" name="height" value="32KP"/>
  <policy domain="resource" name="area" value="900MB"/>
  <policy domain="resource" name="disk" value="4GiB"/>

  This also helped me to significantly speed up the sharpen process !
  
Currently the s/w runs best in WSL (Windows Subsystem for Linux) :-(

It allows you to generate ELF exe files, which you need for memory space and cygwin 
won't allow. Limitation is that it only mounts NTFS HD's and USB devices, NOT a NAS !

It works both in WSL1 and WSL2. WSL2 is much faster in execution time then WSL1.
If you intend to use the tletrack s/w, then once you start a session do:

sudo ntpdate time.microsoft.com 

a few times. This is needed to syncronise the unix clock, which in WSL2 (due to 
hypervirtualisation) can be substantially off from the windows clock.
Currently (Jan 2021) WSL2 and Virtualbox VM's can't run simultaneously. In a 
admin command shell you do:

bcdedit /SET hypervisorlaunchtype off

and reboot, to enable virtualbox VM's. And you do:

bcdedit /SET hypervisorlaunchtype auto

and reboot, to enable WSL2 again.

With some effort in WSL2 it is possible to NFS mount a NAS.

It also runs fine on UBUNTU as an ORACLE Virtualbox VM, but slower.

Install:
- In windows 10 search for "feature"  or "turn windows features" and start "turn windows features on or off"
- Scroll to Windows Subsystem for Linux and enable it, then OK etc etc
- Go to : https://www.microsoft.com/nl-nl/p/ubuntu-1804-lts/9n9tngvndl3q?rtc=1&activetab=pivot:overviewtab
- click download and login with microsoft account - if you don't have it, create one
- You get a number of questions (user, password etc) and end up with a windows command 
  shell in which you type ubuntu1804 and ubuntu starts
- In this shell:
  - sudo apt-get update
  - sudo apt-get install gfortran
  - sudo apt-get install imagemagick
- git clone <Insert the github link here>
- cd weathersat
- Download and unback the resource file from my google drive (see above)
- READ !!!!! the README file
- Define the MMODEL environment variable BEFORE compiling
- chmod 744 buildhrpt.sh
- ./buildhrpt.sh
- ./hrpt.exe --help

Go to the top end try setting the environment variables and run an example


Basic instructions for a virtualbox setup:

- First note that on most computers running a 64 bit VM requires the virtualisation setting 
  in the BIOS to be enabled !
- Install the latest Oracle Virtualbox
- Download the latest Ubuntu LTS ISO with the "amd64" in the title (also supports Intel)
- In virtual box click new and follow all the steps 
  - I would recommen 8 GByte or better for the memory setting
- Finish the VM installation to the end and it will do a reboot
- Start a terminal
- sudo apt-get install gfortran
- sudo apt-get install imagemagick
- sudo apt-get install git
- git clone - get the link from the github page (pulldown button)
- cd weathersat
- Download the google drive file referred to in the README and install here
- uname -a and ensure it has x86_64 as its machine model
- Define the environment variables
- chmod 744 ./buildhrpt.sh
- ./buildhrpt.sh
- and it should work


Enjoy.

Fred Jansen

@redplanet00  (Twitter)

Fred.RocketScientist@gmail.com

====================================================================================
Known Problems section
====================================================================================
-- Corrupted NOAA file (gave no output and also ReadHRPT from David Taylor complained)

I've seen a case of a corrupted NOAA .raw16 file coming from XHRPT. This was the 
recipe for the cure:

1. Dump the file to a log file on a partition visible to windows:

od -A d -t x1 --width=70 NOAA-19_HRPT_20190309_043441.raw16 >> ../Win7/odraw16.log

2. Open the file in wordpad and do a ctrl-F for:

84 02 6F 01 5C 03 9D 01 0F 02 95

if not found, redo step 1 with e.g. a width of 60 or 71 or something

3. Locate the position of the byte before (!) the start of the signature. Let call this 
   number N. Then calculate:
   
   mod(N,22180)    and letś call the result L
   
4. Cut of the dummy info (example) - the + before the number L is crucial :

tail -c +L ../Win7/SDR/To-Process/NOAA-19_HRPT_20190309_043441.raw16 >> NOAA-19_HRPT_20190309_043441.raw16

so in this case L was 3017:

tail -c +3017 ../Win7/SDR/To-Process/NOAA-19_HRPT_20190309_043441.raw16 >> NOAA-19_HRPT_20190309_043441.raw16

5. then run:

./hrpt.exe ./NOAA-19_HRPT_20190309_043441.raw16 S 445 gapcor project=445 gamma=0.9 border debug >& ../Win7/hrpt.log

and check in the log that the records start with the signature hex chain above.
If so, you were succesful

or use a HEX editor, search for the signature and cut out the part before.