# ------------------------------------------------------------
# TERCOR.sml
# ------------------------------------------------------------
# SET WARNING LEVEL: Refer to B5 & B6.
$warnings 3
# ------------------------------------------------------------
# DEFINE VARIABLES RELATED TO TERCOR MODEL:
numeric p,flr;
numeric s,fCB,fBL,fGL,fYL,fRL,fRE;
numeric fNA,fNB,fMA,fMB,fMC,fMD,fME,fMF,fMG;
numeric srfi1,srfi,srfimax;
flr = 0.12;
p = 2.0;
srfimax = 10000;
# ------------------------------------------------------------
# DEFINE VARIABLES FOR GENERAL CHARACTER STRINGS:
# Refer to B7.
string t$,p$,p1$,p2$,p3$,p4$,p5$,p6$,p7$,p8$,p9$;
string p10$,p11$,p12$,p13$,p14$,p15$,p16$,p17$,p18$,p19$;
string p20$,p21$;
# ------------------------------------------------------------
# CLEAR CONSOLE WINDOW & REQUEST REPOSITIONING:
# Refer to B7, B8 & B9.
clear();
p1$ = "CONSOLE-WINDOW ADJUSTMENT\n";
p2$ = "* REPOSITION the CONSOLE WINDOW.\n";
p3$ = "* Then, CLICK the OK Button.";
p$ = p1$ + p2$ + p3$;
PopupMessage(p$);
# ------------------------------------------------------------
# DEFINE proc checkHisto(Band)
# PURPOSE: Check that Input Raster has a Full (Unsampled)
# Histogram. If Not, Then Recompute Histogram.
proc checkHisto (raster Band) begin
local numeric sampleInt;
sampleInt = HistogramGetSampleInterval(Band);
if (sampleInt > 1) begin
DeleteHistogram(Band); CreateHistogram(Band,0);
end
else if (sampleInt == -1) then begin
CreateHistogram(Band,0);
end
end
# ------------------------------------------------------------
# DEFINE PROCEDURE writeTitle: Refer to B13 & B14.
# PURPOSE: WRITES TITLE & AUTHOR INFO TO CONSOLE WINDOW.
proc writeTitle() begin
printf("TERCOR.sml:\n");
printf(" VERSION: November 15, 2005\n");
printf(" PURPOSE: CORRECT SRFI RASTERS ");
printf("FOR TERRAIN-SHADING EFFECTS\n");
printf(" DETAILS: FAQs_by_Jack A & D\n");
printf(" AUTHOR: Dr. Jack F. Paris\n");
printf(" CONTACT INFO: [email protected] ");
printf(" 303-775-1195\n");
printf(" ALLOWED USE: ONLY NON-COMMERCIAL\n\n");
end
# ------------------------------------------------------------
# WRITE TITLE & AUTHOR INFORMATION: Refer to B14.
writeTitle();
# ------------------------------------------------------------
# DECLARE VARIABLES RELATED TO USER INPUTS:
# Refer to B15 & B16.
string site$ = "Stockton, CA";
numeric imager,vnir;
# ------------------------------------------------------------
# DECLARE VARIABLES DIRECTLY RELATED TO USER PROVIDED VALUES:
# Refer to B15 & B16.
string imager$;
# ------------------------------------------------------------
# DECLARE STRING VARIABLE RELATED TO SRFI RASTERS:
# Refer to DXX.
string rtype$ = "16-bit unsigned";
string stype$ = "8-bit unsigned";
# ------------------------------------------------------------
# DECLARE STRING VARIABLE RELATED TO ALL RASTER TYPES:
# Refer to DXX.
# ------------------------------------------------------------
# DECLARE VARIABLES RELATED TO BOOLEAN ENABLERS:
# Refer to B15 & B16
numeric pCB,pBL,pGL,pYL,pRL,pRE;
numeric pNA,pNB,pMA,pMB,pMC,pMD,pME,pMF,pMG,iip,isp,ipp;
# ------------------------------------------------------------
# DECLARE VARIBLES RELATED TO FOR-EACH LOOPS:
numeric i,nlins,ncols,maxi;
# ------------------------------------------------------------
# DECLARE VARIABLES RELATED TO PVI & PBI: Refer to B15 & B16.
numeric pfac,angrad,sinang,cosang;
numeric pbioff,pvioff,srfiNRint,bslineslope,maxpvipbi;
numeric srfiRL,srfiNR,tsrfiNR,pbif,pbi,pvif,pvi;
# ------------------------------------------------------------
# DECLARE LIST OF POSSIBLE INPUT RASTERS:
raster SRFICB1,SRFIBL1,SRFIGL1,SRFIYL1,SRFIRL1,SRFIRE1;
raster SRFINA1,SRFINB1,SRFIMA1,SRFIMB1,SRFIMC1,SRFIMD1;
raster SRFIME1,SRFIMF1,SRFIMG1;
raster PBI1,PVI1,SHADING;
# ------------------------------------------------------------
# DECLARE LIST OF POSSIBLE OUTPUT RASTERS:
raster SRFICB,SRFIBL,SRFIGL,SRFIYL,SRFIRL,SRFIRE;
raster SRFINA,SRFINB,SRFIMA,SRFIMB,SRFIMC,SRFIMD,SRFIME;
raster SRFIMF,SRFIMG;
raster PBI,PVI;
# ------------------------------------------------------------
# DECLARE VARIABLES RELATED TO WAVELENGTH-BASED MODELS:
# Refer to B15 & B16.
numeric wLenCB,wLenBL,wLenGL,wLenYL,wLenRL,wLenRE;
numeric wLenNA,wLenNB,wLenMA,wLenMB,wLenMC,wLenMD,wLenME;
numeric wLenMF,wLenMG;
# ------------------------------------------------------------
# DECLARE VARIABLES RELATED TO SUN ANGLES:
numeric sunelevang,sunaziang;
# ------------------------------------------------------------
# GET NAME OF SITE FROM THE USER: Refer to B17.
t$ = "SITE NAME";
p1$ = "SITE-NAME ENTRY\n";
p2$ = "* ENTER the SITE NAME.\n";
p3$ = "* Then, CLICK the OK Button.\n\n";
p4$ = "SITE-NAME ENTERED:";
p$ = p1$ + p2$ + p3$ + p4$;
site$ = PopupString(p$,site$,t$);
printf(" SITE NAME: %s\n",site$);
# ------------------------------------------------------------
# GET IMAGER NUMBER FROM THE USER: Refer to B17.
p1$ = "IMAGER-NUMBER:\n";
p2$ = " IMAGER\n";
p3$ = " NUMBER: SYSTEM NAME\n";
p4$ = " 1: QuickBird 2 MS\n";
p5$ = " 2: Ikonos 2 MS\n";
p6$ = " 3: OrbView 3 MS\n";
p7$ = " 4: Landsat 7 ETM+\n";
p8$ = " 5: Landsat 5 TM\n";
p9$ = " 6: Landsat 5 MSS\n";
p10$ = " 7: Landsat 4 TM\n";
p11$ = " 8: Landsat 4 MSS\n";
p12$ = " 9: Landsat 3 MSS\n";
p13$ = " 10: Landsat 2 MSS\n";
p14$ = " 11: Landsat 1 MSS\n";
p15$ = " 12: Terra ASTER\n";
p16$ = " 13: Terra MODIS\n";
p17$ = " 14: Aqua MODIS\n";
p18$ = "* Either ACCEPT the Default NUMBER,\n";
p19$ = "* Or, SELECT a Different NUMBER.\n";
p20$ = "* Then, CLICK the OK Button.\n\n";
p21$ = "NUMBER ENTERED:";
p$ = p1$ + p2$ + p3$ + p4$ + p5$ + p6$ + p7$ + p8$ + p9$;
p$ = p$ + p10$ + p11$ + p12$ + p13$ + p14$ + p15$ + p16$;
p$ = p$ + p17$ + p18$ + p19$ + p20$ + p21$;
imager = PopupNum(p$,4,1,14,0);
# ------------------------------------------------------------
# GENERATE IMAGER-SPECIFIC PARAMETERS: Refer to B25.
if (imager == 1) then begin
imager$ = "QuickBird 2 MS";
pBL=1; pNA=1;
wLenBL=0.482; wLenGL=0.548; wLenRL=0.654; wLenNA=0.809;
end
if (imager == 2) then begin
imager$ = "Ikonos 2 MS"; pBL=1; pNA=1;
wLenBL=0.480; wLenGL=0.551; wLenRL=0.665; wLenNA=0.805;
end
if (imager == 3) then begin
imager$ = "OrbView 3 MS"; pBL=1; pNA=1;
wLenBL=0.482; wLenGL=0.548; wLenRL=0.654; wLenNA=0.809;
end
if (imager == 4) then begin
imager$ = "Landsat 7 ETM+"; pBL=1; pNA=1;
pMB=1; pMC=1;
wLenBL=0.482; wLenGL=0.565; wLenRL=0.660; wLenNA=0.825;
wLenMB=1.650; wLenMC=2.220;
end
if (imager == 5) then begin
imager$ = "Landsat 5 TM";
pBL=1; pNA=1; pMB=1; pMC=1;
wLenBL=0.482; wLenGL=0.565; wLenRL=0.660; wLenNA=0.825;
wLenMB=1.650; wLenMC=2.220;
end
if (imager == 6) then begin
imager$ = "Landsat 5 MSS"; pRE=1; pNB=1;
wLenGL=0.550; wLenRL=0.650; wLenRE=0.750;
wLenNB=0.950;
end
if (imager == 7) then begin
imager$ = "Landsat 4 TM";
pBL=1; pNA=1; pMB=1; pMC=1;
wLenBL=0.482; wLenGL=0.548; wLenRL=0.660; wLenNA=0.825;
wLenMB=1.650; wLenMC=2.220;
end
if (imager == 8) then begin
imager$ = "Landast 4 MSS"; pRE=1; pNB=1;
wLenGL=0.550; wLenRL=0.650; wLenRE=0.750;
wLenNB=0.950;
end
if (imager == 9) then begin
imager$ = "Landsat 3 MSS"; pRE=1; pNB=1;
wLenGL=0.550; wLenRL=0.650; wLenRE=0.750;
wLenNB=0.950;
end
if (imager == 10) then begin
imager$ = "Landsat 2 MSS"; pRE=1; pNB=1;
wLenGL=0.550; wLenRL=0.650; wLenRE=0.750;
wLenNB=0.950;
end
if (imager == 11) then begin
imager$ = "Landsat 1 MSS"; pRE=1; pNB=1;
wLenGL=0.550; wLenRL=0.650; wLenRE=0.750;
wLenNB=0.950;
end
if (imager == 12) then begin
imager$ = "Terra ASTER"; pNA=1; pMB=1; pMC=1;
pMD=1; pME=1; pMF=1; pMG=1;
p1$ = "VNIR or ALL OPTION:\n";
p2$ = " OPTION 1: VNIR ONLY\n";
p3$ = " OPTION 2: ALL 9 BANDS\n";
p4$ = "OPTION ENTERED:";
p$ = p1$ + p2$ + p3$ + p4$;
vnir = PopupNum(p$,2,1,2,0);
if (vnir == 1) then begin
pMB=0; pMC=0; pMD=0; pME=0; pMF=0; pMG=0;
end
wLenGL=0.560; wLenRL=0.660; wLenNA=0.820;
wLenMB=1.650; wLenMC=2.165; wLenMD=2.205; wLenME=2.260;
wLenMF=2.330; wLenMG=2.395;
end
if (imager == 13) then begin
imager$ = "Terra MODIS";
pBL=1; pNA=1; pMA=1; pMB=1; pMC=1;
wLenBL=0.469; wLenGL=0.555; wLenRL=0.645; wLenNA=0.8585;
wLenMA=1.240; wLenMB=1.640; wLenMC=2.013;
end
if (imager == 14) then begin
imager$ = "Aqua MODIS";
pBL=1; pNA=1; pMA=1; pMB=1; pMC=1;
wLenBL=0.469; wLenGL=0.555; wLenRL=0.645; wLenNA=0.8585;
wLenMA=1.240; wLenMB=1.640; wLenMC=2.013;
end
printf(" IMAGER: %s\n\n",imager$);
# ------------------------------------------------------------
# GET fGL VALUE FROM USER: Refer to D2.
p1$ = "GL SKY SPECTRAL-IRRADIANCE FRACTION:\n";
p2$ = " fGL ROLE: Determines f for Other Bands)\n";
p3$ = " fGL RANGE: 0.20 to 0.90\n";
p4$ = "* Either ACCEPT the Default FRACTION,\n";
p5$ = "* Or, ENTER a Different FRACTION.\n";
p6$ = "* Then, CLICK the OK Button.\n\n";
p7$ = "FRACTION ENTERED:";
p$ = p1$ + p2$ + p3$ + p4$ + p5$ + p6$ + p7$;
fGL = PopupNum(p$,0.40,0.10,1.0,2);
# ------------------------------------------------------------
# ASSIGN VALUES TO f-FACTORS: Refer to D2.
printf("SKY SPECTRAL-IRRADIANCE FRACTIONS, fBAND:\n");
printf(" BAND wLen fBAND\n");
if (pCB) then begin
fCB = flr + (fGL - flr) * (wLenGL / wLenCB)^p;
printf(" CB%7.3f%7.3f (MODEL)\n",wLenCB,fCB);
end
if (pBL) then begin
fBL = flr + (fGL - flr) * (wLenGL / wLenBL)^p;
printf(" BL%7.3f%7.3f (MODEL)\n",wLenBL,fBL);
end
printf(" GL%7.3f%7.3f (INPUT)\n",wLenGL,fGL);
if (pYL) then begin
fYL = flr + (fGL - flr) * (wLenGL / wLenYL)^p;
printf(" YL%7.3f%7.3f (MODEL)\n",wLenYL,fYL);
end
fRL = flr + (fGL - flr) * (wLenGL / wLenRL)^p;
printf(" RL%7.3f%7.3f (MODEL)\n",wLenRL,fRL);
if (pRE) then begin
fRE = flr + (fGL - flr) * (wLenGL / wLenRE)^p;
printf(" RE%7.3f%7.3f (MODEL)\n",wLenRE,fRE);
end
if (pNA) then begin
fNA = flr + (fGL - flr) * (wLenGL / wLenNA)^p;
printf(" NA%7.3f%7.3f (MODEL)\n",wLenNA,fNA);
end
if (pNB) then begin
fNB = flr + (fGL - flr) * (wLenGL / wLenNB)^p;
printf(" NB%7.3f%7.3f (MODEL)\n",wLenNB,fNB);
end
if (pMA) then begin
fMA = flr + (fGL - flr) * (wLenGL / wLenMA)^p;
printf(" MA%7.3f%7.3f (MODEL)\n",wLenMA,fMA);
end
if (pMB) then begin
fMB = flr + (fGL - flr) * (wLenGL / wLenMB)^p;
printf(" MB%7.3f%7.3f (MODEL)\n",wLenMB,fMB);
end
if (pMC) then begin
fMC = flr + (fGL - flr) * (wLenGL / wLenMC)^p;
printf(" MC%7.3f%7.3f (MODEL)\n",wLenMC,fMC);
end
if (pMD) then begin
fMD = flr + (fGL - flr) * (wLenGL / wLenMD)^p;
printf(" MD%7.3f%7.3f (MODEL)\n",wLenMD,fMD);
end
if (pME) then begin
fME = flr + (fGL - flr) * (wLenGL / wLenME)^p;
printf(" ME%7.3f%7.3f (MODEL)\n",wLenME,fME);
end
if (pMF) then begin
fMF = flr + (fGL - flr) * (wLenGL / wLenMF)^p;
printf(" MF%7.3f%7.3f (MODEL)\n",wLenMC,fMC);
end
if (pMG) then begin
fMG = flr + (fGL - flr) * (wLenGL / wLenMG)^p;
printf(" MC%7.3f%7.3f (MODEL)\n",wLenMC,fMC);
end
printf("\n fBAND MODEL EXPONENT: %6.2f\n",p);
printf(" fBAND MODEL FLOOR: %6.2f\n",flr);
# ------------------------------------------------------------
# GET SUN ELEVATION ANGLE FROM THE USER: Refer to B17.
p1$ = "SUN-ELEVATION-ANGLE ENTRY\n";
p2$ = " UNITS: degrees (above the horizon)\n";
p3$ = " FORMAT: NN.NN\n";
p4$ = " RANGE: 10.00 to 90.00 degrees\n";
p5$ = "* Either ACCEPT the Default ANGLE,\n";
p6$ = "* Or, ENTER a Different ANGLE.\n";
p7$ = "* Then, CLICK the OK Button.\n\n";
p8$ = "SUN-ELEVATION-ANGLE (degrees) ENTERED:";
p$ = p1$ + p2$ + p3$ + p4$ + p5$ + p6$ + p7$ + p8$;
sunelevang = PopupNum(p$,45.18,10,90,2);
printf(" SUN ELEV. ANGLE: %6.2f deg\n",sunelevang);
# ------------------------------------------------------------
# GET SUN AZIMUTH ANGLE FROM THE USER: Refer to EXX.
p1$ = "SUN-AZIMUTH-ANGLE:\n";
p2$ = " UNITS: degrees (relative to true north)\n";
p3$ = " FORMAT: NNN.NN\n";
p4$ = " RANGE: 0.00 to 360.00 degrees\n";
p5$ = "* Either ACCEPT the Default ANGLE,\n";
p6$ = "* Or, ENTER a Different ANGLE.\n";
p7$ = "* Then, CLICK the OK Button.\n\n";
p8$ = "ANGLE ENTERED:";
p$ = p1$ + p2$ + p3$ + p4$ + p5$ + p6$ + p7$ + p8$;
sunaziang = PopupNum(p$,149.41,0,360,2);
printf(" SUN AZI. ANGLE: %6.2f deg\n\n",sunaziang);
# ------------------------------------------------------------
# OPEN INPUT RASTERS: Refer to EXX.
printf("INPUT RASTERS:\n");
if (pCB) then begin
printf(" SRFICB");
GetInputRaster(SRFICB1);
nlins = NumLins(SRFICB1); ncols = NumCols(SRFICB1);
checkHisto(SRFICB1);
end
if (pCB) then begin
if (pBL) then begin
printf(" SRFIBL");
GetInputRaster(SRFIRL1,nlins,ncols,rtype$);
checkHisto(SRFIRL1);
end
end
else begin
if (pBL) then begin
printf(" SRFIBL");
GetInputRaster(SRFIBL1);
nlins = NumLins(SRFIBL1); ncols = NumCols(SRFIBL1);
checkHisto(SRFIBL1);
end
end
printf(" SRFIGL");
if (pBL) then begin
GetInputRaster(SRFIGL1,nlins,ncols,rtype$);
checkHisto(SRFIGL1);
end
else begin
GetInputRaster(SRFIGL1);
nlins = NumLins(SRFIGL1); ncols = NumCols(SRFIGL1);
checkHisto(SRFIGL1);
end
if (pYL) then begin
printf(" SRFIYL");
GetInputRaster(SRFIYL1,nlins,ncols,rtype$);
checkHisto(SRFIYL1);
end
printf(" SRFIRL");
GetInputRaster(SRFIRL1,nlins,ncols,rtype$);
checkHisto(SRFIRL1);
if (pRE) then begin
printf(" SRFIRE");
GetInputRaster(SRFIRE1,nlins,ncols,rtype$);
checkHisto(SRFIRE1);
end
if (pNA) then begin
printf(" SRFINA");
GetInputRaster(SRFINA1,nlins,ncols,rtype$);
checkHisto(SRFINA1);
end
if (pNB) then begin
printf(" SRFINB");
GetInputRaster(SRFINB1,nlins,ncols,rtype$);
checkHisto(SRFINB1);
end
if (pMA) then begin
printf(" SRFIMA");
GetInputRaster(SRFIMA1,nlins,ncols,rtype$);
checkHisto(SRFIMA1);
end
if (pMB) then begin
printf(" SRFIMB");
GetInputRaster(SRFIMB1,nlins,ncols,rtype$);
checkHisto(SRFIMB1);
end
if (pMC) then begin
printf(" SRFIMC");
GetInputRaster(SRFIMC1,nlins,ncols,rtype$);
checkHisto(SRFIMC1);
end
if (pMD) then begin
printf(" SRFIMD");
GetInputRaster(SRFIMD1,nlins,ncols,rtype$);
checkHisto(SRFIMD1);
end
if (pME) then begin
printf(" SRFIME");
GetInputRaster(SRFIME1,nlins,ncols,rtype$);
checkHisto(SRFIME1);
end
if (pMF) then begin
printf(" SRFIMF");
GetInputRaster(SRFIMF1,nlins,ncols,rtype$);
checkHisto(SRFIMF1);
end
if (pMG) then begin
printf(" SRFIMG");
GetInputRaster(SRFIMG1,nlins,ncols,rtype$);
checkHisto(SRFIMG1);
end
printf(" PVI");
GetInputRaster(PVI1,nlins,ncols,rtype$);
checkHisto(PVI1);
printf(" PBI");
GetInputRaster(PBI1,nlins,ncols,rtype$);
checkHisto(PBI1);
printf(" SHADING");
GetInputRaster(SHADING,nlins,ncols,stype$);
checkHisto(SHADING);
# ------------------------------------------------------------
# SET UP OUTPUT RASTERS:
printf("\n\nOUTPUT RASTERS:\n");
if (pCB) then begin
printf(" SRFICB");
GetOutputRaster(SRFICB,nlins,ncols,rtype$);
SetNull(SRFICB,0);
CopySubobjects(SRFICB1,SRFICB,"GEOREF");
CopySubobjects(SRFICB1,SRFICB,"CONTAB");
end
if (pBL) then begin
printf(" SRFIBL");
GetOutputRaster(SRFIBL,nlins,ncols,rtype$);
SetNull(SRFIBL,0);
CopySubobjects(SRFIBL1,SRFIBL,"GEOREF");
CopySubobjects(SRFIBL1,SRFIBL,"CONTAB");
end
printf(" SRFIGL");
GetOutputRaster(SRFIGL,nlins,ncols,rtype$);
SetNull(SRFIGL,0);
CopySubobjects(SRFIGL1,SRFIGL,"GEOREF");
CopySubobjects(SRFIGL1,SRFIGL,"CONTAB");
if (pYL) then begin
printf(" SRFIYL");
GetOutputRaster(SRFIYL,nlins,ncols,rtype$);
SetNull(SRFIYL,0);
CopySubobjects(SRFIYL1,SRFIYL,"GEOREF");
CopySubobjects(SRFIYL1,SRFIYL,"CONTAB");
end
printf(" SRFIRL");
GetOutputRaster(SRFIRL,nlins,ncols,rtype$);
SetNull(SRFIRL,0);
CopySubobjects(SRFIRL1,SRFIRL,"GEOREF");
CopySubobjects(SRFIRL1,SRFIRL,"CONTAB");
if (pRE) then begin
printf(" SRFIRE");
GetOutputRaster(SRFIRE,nlins,ncols,rtype$);
SetNull(SRFIRE,0);
CopySubobjects(SRFIRE1,SRFIRE,"GEOREF");
CopySubobjects(SRFIRE1,SRFIRE,"CONTAB");
end
if (pNA) then begin
printf(" SRFINA");
GetOutputRaster(SRFINA,nlins,ncols,rtype$);
SetNull(SRFINA,0);
CopySubobjects(SRFINA1,SRFINA,"GEOREF");
CopySubobjects(SRFINA1,SRFINA,"CONTAB");
end
if (pNB) then begin
printf(" SRFINB");
GetOutputRaster(SRFINB,nlins,ncols,rtype$);
SetNull(SRFINB,0);
CopySubobjects(SRFINB1,SRFINB,"GEOREF");
CopySubobjects(SRFINB1,SRFINB,"CONTAB");
end
if (pMA) then begin
printf(" SRFIMA");
GetOutputRaster(SRFIMA,nlins,ncols,rtype$);
SetNull(SRFIMA,0);
CopySubobjects(SRFIMA1,SRFIMA,"GEOREF");
CopySubobjects(SRFIMA1,SRFIMA,"CONTAB");
end
if (pMB) then begin
printf(" SRFIMB");
GetOutputRaster(SRFIMB,nlins,ncols,rtype$);
SetNull(SRFIMB,0);
CopySubobjects(SRFIMB1,SRFIMB,"GEOREF");
CopySubobjects(SRFIMB1,SRFIMB,"CONTAB");
end
if (pMC) then begin
printf(" SRFIMC");
GetOutputRaster(SRFIMC,nlins,ncols,rtype$);
SetNull(SRFIMC,0);
CopySubobjects(SRFIMC1,SRFIMC,"GEOREF");
CopySubobjects(SRFIMC1,SRFIMC,"CONTAB");
end
if (pMD) then begin
printf(" SRFIMD");
GetOutputRaster(SRFIMD,nlins,ncols,rtype$);
SetNull(SRFIMD,0);
CopySubobjects(SRFIMD1,SRFIMD,"GEOREF");
CopySubobjects(SRFIMD1,SRFIMD,"CONTAB");
end
if (pME) then begin
printf(" SRFIME");
GetOutputRaster(SRFIME,nlins,ncols,rtype$);
SetNull(SRFIME,0);
CopySubobjects(SRFIME1,SRFIME,"GEOREF");
CopySubobjects(SRFIME1,SRFIME,"CONTAB");
end
if (pMF) then begin
printf(" SRFIMF");
GetOutputRaster(SRFIMF,nlins,ncols,rtype$);
SetNull(SRFIMF,0);
CopySubobjects(SRFIMF1,SRFIMF,"GEOREF");
CopySubobjects(SRFIMF1,SRFIMF,"CONTAB");
end
if (pMG) then begin
printf(" SRFIMG");
GetOutputRaster(SRFIMG,nlins,ncols,rtype$);
SetNull(SRFIMG,0);
CopySubobjects(SRFIMG1,SRFIMG,"GEOREF");
CopySubobjects(SRFIMG1,SRFIMG,"CONTAB");
end
printf(" PVI");
GetOutputRaster(PVI,nlins,ncols,rtype$);
SetNull(PVI,0);
CopySubobjects(PVI1,PVI,"GEOREF");
CopySubobjects(PVI1,PVI,"CONTAB");
CloseRaster(PVI1);
printf(" PBI");
GetOutputRaster(PBI,nlins,ncols,rtype$);
SetNull(PBI,0);
CopySubobjects(PBI1,PBI,"GEOREF");
CopySubobjects(PBI1,PBI,"CONTAB");
CloseRaster(PBI1);
# ------------------------------------------------------------
# PRODUCE OUTPUT RASTERS: Refer to EXX.
printf("\n\n");
printf("CORRECTING SRFT RASTERS:\n");
array numeric tc[256];
for i=0 to 254 begin
tc[i] = i / 180;
end
# CB DATA:
if (pCB) then begin
printf(" SRFICB");
for each SRFICB1 begin
iip = 1; isp = 1;
if (IsNull(SRFICB1)) then iip = 0;
if (IsNull(SHADING)) then isp = 0;
ipp = iip * isp;
if (ipp) then begin
srfi1 = SRFICB1; s = SHADING; s = tc[s];
srfi = round(srfi1/(s*(1-fCB)+fCB));
if (srfi < 1) then srfi = 1;
if (srfi > srfimax) then srfi = srfimax;
SRFICB = srfi;
end
end
CreateHistogram(SRFICB,0); CreatePyramid(SRFICB,0);
CloseRaster(SRFICB1); CloseRaster(SRFICB);
end
# BL DATA:
if (pBL) then begin
printf(" SRFIBL");
for each SRFIBL1 begin
iip = 1; isp = 1;
if (IsNull(SRFIBL1)) then iip = 0;
if (IsNull(SHADING)) then isp = 0;
ipp = iip * isp;
if (ipp) then begin
srfi1 = SRFIBL1; s = SHADING; s = tc[s];
srfi = round(srfi1/(s*(1-fBL)+fBL));
if (srfi < 1) then srfi = 1;
if (srfi > srfimax) then srfi = srfimax;
SRFIBL = srfi;
end
end
CreateHistogram(SRFIBL,0); CreatePyramid(SRFIBL,0);
CloseRaster(SRFIBL1); CloseRaster(SRFIBL);
end
# GL DATA:
printf(" SRFIGL");
for each SRFIGL1 begin
iip = 1; isp = 1;
if (IsNull(SRFIGL1)) then iip = 0;
if (IsNull(SHADING)) then isp = 0;
ipp = iip * isp;
if (ipp) then begin
srfi1 = SRFIGL1; s = SHADING; s = tc[s];
srfi = round(srfi1/(s*(1-fGL)+fGL));
if (srfi < 1) then srfi = 1;
if (srfi > srfimax) then srfi = srfimax;
SRFIGL = srfi;
end
end
CreateHistogram(SRFIGL,0); CreatePyramid(SRFIGL,0);
CloseRaster(SRFIGL1); CloseRaster(SRFIGL);
# YL DATA:
if (pYL) then begin
printf(" SRFIYL");
for each SRFIYL1 begin
iip = 1; isp = 1;
if (IsNull(SRFIYL1)) then iip = 0;
if (IsNull(SHADING)) then isp = 0;
ipp = iip * isp;
if (ipp) then begin
srfi1 = SRFIYL1; s = SHADING; s = tc[s];
srfi = round(srfi1/(s*(1-fYL)+fYL));
if (srfi < 1) then srfi = 1;
if (srfi > srfimax) then srfi = srfimax;
SRFIYL = srfi;
end
end
CreateHistogram(SRFIYL,0); CreatePyramid(SRFIYL,0);
CloseRaster(SRFIYL1); CloseRaster(SRFIYL);
end
# RL DATA:
printf(" SRFIRL");
for each SRFIRL1 begin
iip = 1; isp = 1;
if (IsNull(SRFIRL1)) then iip = 0;
if (IsNull(SHADING)) then isp = 0;
ipp = iip * isp;
if (ipp) then begin
srfi1 = SRFIRL1; s = SHADING; s = tc[s];
srfi = round(srfi1/(s*(1-fRL)+fRL));
if (srfi < 1) then srfi = 1;
if (srfi > srfimax) then srfi = srfimax;
SRFIRL = srfi;
end
end
CreateHistogram(SRFIRL,0); CreatePyramid(SRFIRL,0);
CloseRaster(SRFIRL1);
# RE DATA:
if (pRE) then begin
printf(" SRFIRE");
for each SRFIRE1 begin
iip = 1; isp = 1;
if (IsNull(SRFIRE1)) then iip = 0;
if (IsNull(SHADING)) then isp = 0;
ipp = iip * isp;
if (ipp) then begin
srfi1 = SRFIRE1; s = SHADING; s = tc[s];
srfi = round(srfi1/(s*(1-fRE)+fRE));
if (srfi < 1) then srfi = 1;
if (srfi > srfimax) then srfi = srfimax;
SRFIRE = srfi;
end
end
CreateHistogram(SRFIRE,0); CreatePyramid(SRFIRE,0);
CloseRaster(SRFIRE1); CloseRaster(SRFIRE);
end
# NA DATA:
if (pNA) then begin
printf(" SRFINA");
for each SRFINA1 begin
iip = 1; isp = 1;
if (IsNull(SRFINA1)) then iip = 0;
if (IsNull(SHADING)) then isp = 0;
ipp = iip * isp;
if (ipp) then begin
srfi1 = SRFINA1; s = SHADING; s = tc[s];
srfi = round(srfi1/(s*(1-fNA)+fNA));
if (srfi < 1) then srfi = 1;
if (srfi > srfimax) then srfi = srfimax;
SRFINA = srfi;
end
end
CreateHistogram(SRFINA,0); CreatePyramid(SRFINA,0);
CloseRaster(SRFINA1);
end
# NB DATA:
if (pNB) then begin
printf(" SRFINB");
for each SRFINB1 begin
iip = 1; isp = 1;
if (IsNull(SRFINB1)) then iip = 0;
if (IsNull(SHADING)) then isp = 0;
ipp = iip * isp;
if (ipp) then begin
srfi1 = SRFINB1; s = SHADING; s = tc[s];
srfi = round(srfi1/(s*(1-fNB)+fNB));
if (srfi < 1) then srfi = 1;
if (srfi > srfimax) then srfi = srfimax;
SRFINB = srfi;
end
end
CreateHistogram(SRFINB,0); CreatePyramid(SRFINB,0);
CloseRaster(SRFINB1);
end
# MA DATA:
if (pMA) then begin
printf(" SRFIMA");
for each SRFIMA1 begin
iip = 1; isp = 1;
if (IsNull(SRFIMA1)) then iip = 0;
if (IsNull(SHADING)) then isp = 0;
ipp = iip * isp;
if (ipp) then begin
srfi1 = SRFIMA1; s = SHADING; s = tc[s];
srfi = round(srfi1/(s*(1-fMA)+fMA));
if (srfi < 1) then srfi = 1;
if (srfi > srfimax) then srfi = srfimax;
SRFIMA = srfi;
end
end
CreateHistogram(SRFIMA,0); CreatePyramid(SRFIMA,0);
CloseRaster(SRFIMA1); CloseRaster(SRFIMA);
end
# MB DATA:
if (pMB) then begin
printf(" SRFIMB");
for each SRFIMB1 begin
iip = 1; isp = 1;
if (IsNull(SRFIMB1)) then iip = 0;
if (IsNull(SHADING)) then isp = 0;
ipp = iip * isp;
if (ipp) then begin
srfi1 = SRFIMB1; s = SHADING; s = tc[s];
srfi = round(srfi1/(s*(1-fMB)+fMB));
if (srfi < 1) then srfi = 1;
if (srfi > srfimax) then srfi = srfimax;
SRFIMB = srfi;
end
end
CreateHistogram(SRFIMB,0); CreatePyramid(SRFIMB,0);
CloseRaster(SRFIMB1); CloseRaster(SRFIMB);
end
# MC DATA:
if (pMC) then begin
printf(" SRFIMC");
for each SRFIMC1 begin
iip = 1; isp = 1;
if (IsNull(SRFIMC1)) then iip = 0;
if (IsNull(SHADING)) then isp = 0;
ipp = iip * isp;
if (ipp) then begin
srfi1 = SRFIMC1; s = SHADING; s = tc[s];
srfi = round(srfi1/(s*(1-fMC)+fMC));
if (srfi < 1) then srfi = 1;
if (srfi > srfimax) then srfi = srfimax;
SRFIMC = srfi;
end
end
CreateHistogram(SRFIMC,0); CreatePyramid(SRFIMC,0);
CloseRaster(SRFIMC1); CloseRaster(SRFIMC);
end
# MD DATA:
if (pMD) then begin
printf(" SRFIMD");
for each SRFIMD1 begin
iip = 1; isp = 1;
if (IsNull(SRFIMD1)) then iip = 0;
if (IsNull(SHADING)) then isp = 0;
ipp = iip * isp;
if (ipp) then begin
srfi1 = SRFIMD1; s = SHADING; s = tc[s];
srfi = round(srfi1/(s*(1-fMD)+fMD));
if (srfi < 1) then srfi = 1;
if (srfi > srfimax) then srfi = srfimax;
SRFIMD = srfi;
end
end
CreateHistogram(SRFIMD,0); CreatePyramid(SRFIMD,0);
CloseRaster(SRFIMD1); CloseRaster(SRFIMD);
end
# ME DATA:
if (pME) then begin
printf(" SRFIME");
for each SRFIME1 begin
iip = 1; isp = 1;
if (IsNull(SRFIME1)) then iip = 0;
if (IsNull(SHADING)) then isp = 0;
ipp = iip * isp;
if (ipp) then begin
srfi1 = SRFIME1; s = SHADING; s = tc[s];
srfi = round(srfi1/(s*(1-fME)+fME));
if (srfi < 1) then srfi = 1;
if (srfi > srfimax) then srfi = srfimax;
SRFIME = srfi;
end
end
CreateHistogram(SRFIME,0); CreatePyramid(SRFIME,0);
CloseRaster(SRFIME1); CloseRaster(SRFIME);
end
# MF DATA:
if (pMF) then begin
printf(" SRFIMF");
for each SRFIMF1 begin
iip = 1; isp = 1;
if (IsNull(SRFIMF1)) then iip = 0;
if (IsNull(SHADING)) then isp = 0;
ipp = iip * isp;
if (ipp) then begin
srfi1 = SRFIMF1; s = SHADING; s = tc[s];
srfi = round(srfi1/(s*(1-fMF)+fMF));
if (srfi < 1) then srfi = 1;
if (srfi > srfimax) then srfi = srfimax;
SRFIMF = srfi;
end
end
CreateHistogram(SRFIMF,0); CreatePyramid(SRFIMF,0);
CloseRaster(SRFIMF1); CloseRaster(SRFIMF);
end
# MG DATA:
if (pMG) then begin
printf(" SRFIMG");
for each SRFIMG1 begin
iip = 1; isp = 1;
if (IsNull(SRFIMG1)) then iip = 0;
if (IsNull(SHADING)) then isp = 0;
ipp = iip * isp;
if (ipp) then begin
srfi1 = SRFIMG1; s = SHADING; s = tc[s];
srfi = round(srfi1/(s*(1-fMG)+fMG));
if (srfi < 1) then srfi = 1;
if (srfi > srfimax) then srfi = srfimax;
SRFIMG = srfi;
end
end
CreateHistogram(SRFIMG,0); CreatePyramid(SRFIMG,0);
CloseRaster(SRFIMG1); CloseRaster(SRFIMG);
end
CloseRaster(SHADING);
printf("\n\n");
# ------------------------------------------------------------
# PRODUCE VALUES FOR PVI & PBI RASTERS:
printf("PRODUCING PVI & PBI RASTERS.");
srfiNRint = 254; bslineslope = 1.086;
pbioff = 100; pvioff = 1000;
pfac = 0.2723659; maxpvipbi = 3000;
angrad = -atan(bslineslope);
sinang = sin(angrad); cosang = cos(angrad);
for each SRFIRL begin
srfiRL = SRFIRL;
if (srfiRL > 0) then begin
if (pNA) then srfiNR = SRFINA;
if (pNB and pNA==0) then begin
srfiNR = SRFINB;
end
tsrfiNR = srfiNR - srfiNRint;
pbif = srfiRL * cosang - srfiNR * sinang;
pbi = round(pfac * pbif)+ pbioff;
pvif = srfiRL * sinang + srfiNR * cosang;
pvi = round(pfac * pvif) + pvioff;
if (pbi < 1) then pbi = 1;
if (pbi > maxpvipbi) then pbi = maxpvipbi;
if (pvi < 1) then pvi = 1;
if (pvi > maxpvipbi) then pvi = maxpvipbi;
PBI = pbi; PVI = pvi;
end
end
CreateHistogram(PBI,0); CreateHistogram(PVI,0);
CreatePyramid(PBI,0); CreatePyramid(PVI,0);
CloseRaster(PBI); CloseRaster(PVI);
if (pNA) then CloseRaster(SRFINA);
if (pNB) then CloseRaster(SRFINB);
CloseRaster(SRFIRL);
# ------------------------------------------------------------
printf("TO SAVE THE CONSOLE WINDOW TEXT AS A REPORT:\n");
printf(" 1. RIGHT-CLICK IN THE CONSOLE WINDOW.\n");
printf(" 2. SELECT THE Save As... OPTION.\n");
printf(" 3. NAVIGATE TO THE DESIRED LOCATION.\n");
printf(" 4. PROVIDE A REPORT NAME (or OVERWRITE).\n");
printf(" 5. CLICK OK.");
Download