= Overview = Archiving SPC-focus products. There are two major categories (global and regional) with two subdivisions (topography and backplanes -- aka other data). The regional products are derived from bigmaps. === Tools === These tools are stored in GitHub as [[https://github.com/StereoPhotoClinometry/archivingTools|archivingTools]]. * [[convertLatLon]] - Reads an X,Y,Z (or X,Y,Z,V) list and converts it into lat, lon, distance or lat, lon, value * [[backplanesGMT]] - Reads a bigmap and generates a list of latlon, radius, albedo, slope. Use paste latlon radius > tmp for GMT. * [[backplanesISIS]] - Reads a bigmap and generates 2D matrix for radius, albedo, slope, lat, lon * [[phasei]] - Reads a bigmap and image (the SUMFILE) and generates 2D matrix for incidence, emission, phase angles * coverage_p_low -- Reads all images and maplets. Creates a pgm that has every pixel marked for every image that it has, counter is only 1. * [[map_coverage_p_low]] - Counts the number of images that cover each bigmap pixel. Creates MAPNM-cov.txt and MAPNM-cov.pgm * [[bustGrid]] - Takes an ascii grid output and puts each entry on its own line. Needed to take grid and paste it with lat/lon to make grid -------- = Grid to GeoTIF = {{{ var=radius.ll new=big3 argR="-R280/287/-9/-3" argR="-R0/360/-90/90" argI="-I0.04" date echo "Blockmean" gmt blockmean $var $argR $argI > out echo "sphereinterpolate" cat out | gmt sphinterpolate $argR $argI -Gtest.nc echo "gdal translate" gdal_translate -of GTiff -b 1 -a_srs support/60300.prj NETCDF:test.nc $new.tif }}} ------ = Global = Global ICQ resolution to bin sizes [degrees] || ICQ Q Value || Vertices || Bin Size [deg] || || 512 || 1,579,014 || 0.041 || || 256 || 396,294 || 0.164 || || 128 || 99,846 || 0.649 || || 64 || 25,350 || 2.556 || || 32 || 6,534 || 9.917 || == Topography == {{{ /bin/cp -f SHAPE.TXT tmp vi tmp 1,%s/D/E/g ZZ convertLatLon tmp > radius.ll }}} == Sigmas == {{{ convertLatLon tmp 1 > sigma.ll }}} * Use sigma.ll as the input to GMT == Number of Images == * Run globalStat * Use global-cov.ll as the input to GMT {{{ var=global-cov.ll new=tCov argR="-R0/360/-90/90" argI="-I1" }}} == Best Resolution == * Run globalStat * Use global-res as the input to GMT {{{ var=global-res.ll new=tRes argR="-R0/360/-90/90" argI="-I1" }}} == Slope == == Albedo == {{{ /bin/cp -f SHAPEA.TXT tmp vi tmp 1,$s/D/E/g ZZ convertLatLon tmp 1 > albedo.ll }}} -------- = Regional = Run backplanesGMT (which is in ORExSPCsupport) (or backplanesISIS) == Topography == * Run backplanesGMT {{{ # Make the backplanes name=LEADEQ echo $name | ~/bin/backplanesGMT }}} Output: {{{ Version: 1.50000000 Input map name (only 6 char no MAPFILES and .MAP) BIGMAP position V 120.38036346435547 -510.48150634765625 -54.750602722167969 Ux -2.9186813160777092E-002 0.12163921445608139 -0.99214518070220947 Uy 0.97239929437637329 0.23332308232784271 0.0000000000000000 Uz 0.23149037361145020 -0.96476125717163086 -0.12509183585643768 QSZ 99 Lon 279.612366 286.917664 Lat -9.55961323 -2.35706925 }}} We'll need the Lon and Lat range for "argR", and QSZ will give an idea of what you need for "argI". {{{ # Convert them into lon/lat/val paste $name-lonlat.txt $name-r.txt > tmp.ll var=tmp.ll argR="-R279.612366/286.917664/-9.55961323/-2.35706925" argI="-I199+n" # Run the rebinning/interpolating routines echo "Blockmean" gmt blockmean $var $argR $argI > out echo "sphereinterpolate" cat out | gmt sphinterpolate $argR $argI -Gtest.nc }}} Now you need some math to convert the Lon range number of nodes (i.e. "argI") to deg/px. (286.917664 - 279.612366) / 199 = 0.03671 {{{ # Convert to geoTiff echo "gdal translate" degPx=0.0367 gdal_translate -of ISIS3 -tr $degPx $degPx -r bilinear -b 1 -co TARGET_NAME=Tethys -co DATA_LOCATION=GEOTIFF -a_srs support/60300.prj NETCDF:test.nc $name-r-ISIS.lbl }}} If you don't want to make it readable by ISIS, and don't care about square pixels, you can run the following {{{ gdal_translate -of GTiff -b 1 -a_srs support/60300.prj NETCDF:test.nc $name-r.tif }}} prj files can be found [[https://spatialreference.org/ref/iau2000/|here]] === How to trim size === * Use ArcMap to find the new lat and lon ranges * Remake geotiff with new lat and lon ranges * Measure the size of the lat range. Divide this number by the GSD from SPC to get number of pixels. * Divide number of pixels by lon range to get degPx * Update config.txt with these new values (i.e. lat and lon ranges, number of pixels across, degPx) * Make geotiff again == Slope == * See topography == Albedo == * See topography == Sigmas == * Run backplanesGMT (this gets lon/lat only file) * Take the SIGMAS.pgm and convert it to lon/lat/val files {{{ name=LEADEQ std2isis from=SIGMAS.pgm to=tmp.cub isis2ascii from=tmp.cub to=tmp1.txt sed '1,2d' tmp1.txt > tmp2.txt ~/bin/bustGrid tmp2.txt > $name-sig.txt paste $name-lonlat.txt $name-sig.txt > sigma.ll var=sigma.ll }}} * Continue with topography script == Number of Images == * run map_coverage_p_low * Convert coverage out to lon/lat/value {{{ name=LEADEQ bustGrid $name-cov.txt > tmp1.txt paste $name-lonlat.txt tmp1.txt > imgNum.ll var=imgNum.ll }}} * Continue with topography script == Photometric Data == This includes emission angle, incidence angle and phase angle. Run phasei (located in ORExSPCsupport) Run process similar to Number of Images (i.e. use "bustGrid" and "paste") == Running GMT on the PSI Hawk Cluster == To begin making projected geoTiffs for ArcGIS and ISIS, start with the txt files generated when making the photometric cubes. The scripts used produce the following txt files: *IoverF--.txt (Ex. IoverF-LEADHL-N1563651588.txt) *LEADHL-a.TXT *LEADHL-e.TXT *LEADHL-i.TXT The above files need to be uploaded to the Hawk cluster, while keeping their parent directory structure. Go to the photocube directory for each bigmap {{{ rsync -hapvP --exclude=*.cub --exclude=*.prt N1* jweirich@hawk.psi.edu:/home/jweirich/Tethys/geoTiffTethys/LEADHL/ }}} Also copy the output from backplanesGMT and config.txt Ex. config.txt {{{ BIGMAP=LEADHL ARGR=-R226.209259/261.399750/3.66536975/33.9945526 ARGI=-I801+n DEGPX=0.0378641483770287 PROJECT=60300 BODY=Tethys RUNALL=No VERSION=1 }}} The basic command being run is similar to Topography step above {{{ echo "Blockmean" gmt blockmean $var $argR $argI > ${product}out echo "sphereinterpolate" cat ${product}out | gmt sphinterpolate $argR $argI -G${product}.nc }}} where $product is going to be IoverF, a, e, or i. Everything else is just a wrapper to deconflict all the files since we're running in parallel. I have the GMT commands above in a script called makeaGeoBig.sh which reads most of the inputs from config.txt. makeGeoBig.sh is called by makeIoverFGeoBig.sh, makePhaseGeoBig.sh, makeEmissGeoBig.sh, and makeIncGeoBig.sh The key lines in these scripts are: (Ex. below is for Incidence) {{{ echo "getting ready for GMT ..." ~/bin/bustGrid $bigmap-i.TXT > bustedi-$bigmap-$img.txt paste ../$bigmap-lonlat.txt bustedi-$bigmap-$img.txt > $var echo "calling makeGeoTiffBig.sh" sh ~/bin/makeGeoTiffBig.sh $var $product }}} $product is as before, while $var is IoverF.ll, phase.ll, emiss.ll, or inc.ll The above four scripts are called by run4.sh, which is essentially a copy/paste script with different images. Note the use of & and "wait". & is to get multiple processors running on each machine, and "wait" is to keep the script running until all four processes complete. The key lines from run4.sh look like this: {{{ cd $img sh ~/bin/makeIoverFGeoBig.sh $img & sh ~/bin/makePhaseGeoBig.sh $img & sh ~/bin/makeEmissGeoBig.sh $img & sh ~/bin/makeIncGeoBig.sh $img & cd .. wait }}} Since there are 8 processors per Hawk node, we want to call 2 images per node. To do this we have a script that uses sbatch commands. Here's what one of those looks like: {{{ #!/bin/bash # #SBATCH -J Call1.sh # Job name #SBATCH -t 2-0:00 # time (D-HH:MM) #SBATCH --partition=cpu #SBATCH -o slurm.%N.%j.out # STDOUT #SBATCH -e slurm.%N.%j.err # STDERR sh ~/bin/run4.sh N1563651588 & sh ~/bin/run4.sh N1563651648 & wait }}} Rather than type each by hand, I instead have a script (makeScript.sh) to make a bunch of these: {{{ # 29 Sep 2021 - John R. Weirich # Script to call different photo scripts imgList=$1 if [ -z $imgList ]; then echo "Please select a list of images" echo "Usage: " exit fi list=`cat $imgList` top=1 count=0 rm -f Call*.sh for i in $list do if [ $top = "1" ] then let count=$count+1 out="Call$count.sh" echo "#!/bin/bash" > $out echo "#" >> $out echo "#SBATCH -J $out # Job name" >> $out echo "#SBATCH -t 2-0:00 # time (D-HH:MM)" >> $out echo "#SBATCH --partition=cpu" >> $out echo "#SBATCH -o slurm.%N.%j.out # STDOUT" >> $out echo "#SBATCH -e slurm.%N.%j.err # STDERR" >> $out echo " " >> $out echo "sh ~/bin/run4.sh $i &" >> $out top=2 else echo "sh ~/bin/run4.sh $i &" >> $out echo " " >> $out echo "wait" >> $out top=1 fi done if [ $top = "2" ] then echo " " >> $out echo "wait" >> $out fi }}} To fire everything off, you'll make a copy/paste that will look like this: {{{ sbatch Call1.sh sbatch Call2.sh etc. }}} -------- = Make Pretties = {{{ gmt begin GMT_cont gmt set GMT_THEME cookbook gmt grdcontour test.nc gmt end show #gmt grdcontour test.nc -C10 -A50 gmt begin GMT_img gmt set GMT_THEME cookbook #gmt makecpt -Crainbow gmt grdimage test.nc -JM6i -B -BWSnE gmt colorbar -DJTC -Bxa -By+lm gmt end show gmt begin GMT_img gmt makecpt -Crainbow gmt set GMT_THEME cookbook gmt grdimage test.nc $argR -JM6i -B -BWSnE #gmt colorbar -DJTC -I0.4 -Bxa -By+lm gmt colorbar -DJTC -Bxa -By+lm gmt end show }}} == Setup == Use Dropbox spcShare directory {{{ cd /opt/local/spc/ relink ~/Dropbox/spcShare/Tethys . }}} Then in the working directory {{{ relink.sh /opt/local/spc/Tethys/support . relink.sh /opt/local/spc/Tethys/mapConfig . }}}