# ------------------------------------------------------------
# WATER.sml
# ------------------------------------------------------------
# SET WARNING LEVEL:
$warnings 3
# ------------------------------------------------------------
# DEFINE PROCEDURE writeTitle:
# PURPOSE: WRITES TITLE & AUTHOR INFO TO CONSOLE WINDOW.
proc writeTitle() begin
printf("WATER.sml:\n");
printf(" VERSION: November 15, 2005\n");
printf(" PURPOSE: MAKE A MERGED WATER & ");
printf("LAND IMAGE.\n");
printf(" DETAILS: FAQs_by_Jack G\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
# ------------------------------------------------------------
# DECLARE VARIABLES:
numeric nlins,ncols;
numeric m,d,tLP,v,nc,wm,pW,pL,x,y,yt,facBL,offNA;
numeric sBL_8,sGL_8,sRL_8,sNA_8;
numeric srfiBL,srfiNA;
array sBL_8[65536],sGL_8[65536],sRL_8[65536],sNA_8[65536];
numeric i,mini,srfiBLmaxW,srfiGLmaxW,srfiRLmaxW;
numeric srfiBLmaxL,srfiGLmaxL,srfiRLmaxL,srfiNAmaxL;
numeric redderW,brighterW,brightenDarkW;
numeric redderL,brighterL,brightenDarkL;
numeric pW,pL;
string stype$;
raster SRFIBL,SRFIGL,SRFIRL,SRFINA,WATERMASK;
raster IMAGE;
numeric iB,iG,iR,fB,fG,fR,fCLR;
# ------------------------------------------------------------
# DEFINE VALUES FOR fR, and fG:
fR = 2^16; fG = 2^8; fB = 1;
# ------------------------------------------------------------
# DEFINE STRING VARIABLES FOR POPUP WINDOWS:
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:
p1$ = "CONSOLE-WINDOW ADJUSTMENT:\n";
p2$ = "* REPOSITION the CONSOLE WINDOW.\n";
p3$ = "* Then, CLICK OK.";
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
# ------------------------------------------------------------
# USER INPUTS:
# ENHANCEDMENT DEFAULTS:
redderW = 1.18;
brighterW = 1.00;
brightenDarkW = 1.43;
redderL = 0.95;
brighterL = 1.00;
brightenDarkL = 1.25;
# ------------------------------------------------------------
# TNTmips VERSION / PATCH INPUT:
p1$ = "VERSION / PATCH ENTRY:\n";
p2$ = " Yes: Running TNTmips V7.1 With PATCH\n";
p3$ = " Dated July 13, 2005, or LATER\n";
p4$ = " No: Running Earlier VERSION or PATCH.";
p$ = p1$ + p2$ + p3$ + p4$;
fCLR = PopupYesNo(p$,1);
if (fCLR == 1) then begin
fR = 1; fB = 2^16;
end
# ------------------------------------------------------------
# ACCEPT ENHANCEMENT DEFAULTS:
p$ = "Do You Want to Apply Default ENHANCEMENTS?";
d = PopupYesNo(p$,0);
# ------------------------------------------------------------
# WATER-MASK OPTION:
p1$ = "WATER-MASK OPTION ENTRY:\n";
p2$ = " OPTIONS:\n";
p3$ = " 1: Create a NEW WATERMASK Raster.\n";
p4$ = " 2: Use an EXISTING WATERMASK Raster.\n";
p5$ = " Existing Raster Must Match SRFI Rasters.\n";
p6$ = "OPTION ENTERED:";
p$ = p1$ + p2$ + p3$ + p4$ + p5$ + p6$;
m = PopupNum(p$,2,1,2,0);
# ------------------------------------------------------------
# LAND-IMAGE OPTION:
p1$ = "LAND-IMAGE OPTION ENTRY:\n";
p2$ = " OPTIONS:\n";
p3$ = " 1: COLOR INFRARED IMAGE (CIR).\n";
p4$ = " 2. NATURAL COLOR IMAGE.\n";
p5$ = " 3: GRAYSCALE (GREEN BAND ONLY) IMAGE.\n";
p6$ = " 4: GRAYSCALE (NIR BAND ONLY) IMAGE.\n";
p7$ = "OPTION ENTERED:";
p$ = p1$ + p2$ + p3$ + p4$ + p5$ +p6$ +p7$;
tLP = PopupNum(p$,4,1,4,0);
# ------------------------------------------------------------
# MAKE WATER-COLORS REDDER:
p1$ = "REDDER-WATER-COLORS PARAMETER ENTRY:\n";
p2$ = "VALUE OF redderW ENTERED:";
p$ = p1$ + p2$;
if (d == 0) then begin
redderW = PopupNum(p$,redderW,0.1,2,2);
end
# ------------------------------------------------------------
# MAKE ALL-WATER-COLORS BRIGHTER:
p1$ = "BRIGHTER-WATER-COLORS PARAMETER ENTRY:\n";
p2$ = "VALUE FOR brighterW ENTERED:";
p$ = p1$ + p2$;
if (d == 0) then begin
brighterW = PopupNum(p$,brighterW,0.1,2,2);
end
# ------------------------------------------------------------
# MAKE DARK-WATER-COLORS BRIGHTER:
p1$ = "BRIGHTEN-DARK-WATER-COLORS PARAMETER ENTRY:\n";
p2$ = "VALUE OF brightenDarkW ENTERED:";
p$ = p1$ + p2$;
if (d == 0) then begin
brightenDarkW = PopupNum(p$,brightenDarkW,0.1,2,2);
end
# ------------------------------------------------------------
# MAKE LAND-COLORS REDDER:
p1$ = "REDDER-LAND-COLORS PARAMETER ENTRY:\n";
p2$ = "VALUE OF redderL ENTERED:";
p$ = p1$ + p2$;
if (d == 0) then begin
if (tLP < 3) then begin
redderL = PopupNum(p$,redderL,0.1,2,2);
end
end
# ------------------------------------------------------------
# MAKE ALL-LAND-COLORS BRIGHTER:
p1$ = "BRIGHTER-LAND-COLORS PARAMETER ENTRY:\n";
p2$ = "VALUE OF brighterL ENTERED:";
p$ = p1$ + p2$;
if (d == 0) then begin
brighterL = PopupNum(p$,brighterL,0.1,2,2);
end
# ------------------------------------------------------------
# MAKE DARK-LAND-COLORS BRIGHTER:
p1$ = "BRIGHTEN-DARK-LAND-COLORS PARAMETER ENTRY:\n";
p2$ = "VALUE OF brightenDarkL ENTERED:";
p$ = p1$ + p2$;
if (d == 0) then begin
brightenDarkL = PopupNum(p$,brightenDarkL,0.1,2,2);
end
# ------------------------------------------------------------
# LAND-WATER BOUNDARY CONTROL PARAMETERS:
if (m == 1) then begin
p1$ = "LINE-EQUATION MULTIPLIER ENTRY:\n";
p2$ = " RANGE: 0.00 to 10.00\n";
p3$ = "VALUE OF facBL ENTERED:";
p$ = p1$ + p2$ + p3$;
facBL = PopupNum(p$,0.16,0.0,10,2);
p1$ = "LINE-EQUATION SRFINA OFFSET ENTRY:\n";
p2$ = " RANGE: 1 to 2000\n";
p3$ = "VALUE OF offNA ENTERED:";
p$ = p1$ + p2$ + p3$;
offNA = PopupNum(p$,925,1,3000,0);
end
# ------------------------------------------------------------
# WRITE TITLE & AUTHOR INFORMATION:
clear();
writeTitle();
# ------------------------------------------------------------
# PRINT USER INPUTS TO CONSOLE WINDOW:
if (m == 1) then begin
printf("A NEW WATER-MASK IS BEING CREATED.\n\n");
end
if (m == 2) then begin
printf("AN EXISTING WATER-MASK IS BEING USED.\n");
printf("NOTE: Existing WATER-MASK Must Match");
printf(" SRFI rasters.\n\n");
end
if (tLP == 1) then begin
printf("LAND PIXELS WILL HAVE CIR COLORS.");
end
if (tLP == 2) then begin
printf("LAND PIXELS WILL HAVE NATURAL COLORS.");
end
if (tLP == 3) then begin
printf("LAND PIXELS WILL HAVE GREEN-BAND GRAY TONES.");
end
if (tLP == 4) then begin
printf("LAND PIXELS WILL HAVE NIR-BAND GRAY TONES.");
end
printf("\n\n");
printf("IMAGE ENHANCEMENT PARAMETERS FOR WATER PIXELS:\n");
printf(" redderW = %4.2f\n",redderW);
printf(" brighterW = %4.2f\n",brighterW);
printf(" brightenDarkW = %4.2f\n\n",brightenDarkW);
printf("IMAGE ENHANCEMENT PARAMETERS FOR LAND PIXELS:\n");
if (tLP < 3) then begin
printf(" redderL = %4.2f\n",redderL);
end
printf(" brighterL = %4.2f\n",brighterL);
printf(" brightenDarkL = %4.2f\n",brightenDarkL);
printf("\n");
if (m == 1) then begin
printf("WATERMASK CONTROL PARAMETERS:\n");
printf(" facBL = %4.2f\n",facBL);
printf(" offNA = %4d\n\n",offNA);
end
# ------------------------------------------------------------
# OPEN SRFIsfc RASTERS:
printf("OPEN SRFIsfc RASTERS (ANY APPROPRIATE SOURCE):\n");
printf(" WARNING! MUST be SURFACE-RELATED.\n");
printf(" SRFIBL");
GetInputRaster(SRFIBL);
nlins = NumLins(SRFIBL); ncols = NumCols(SRFIBL);
stype$ = RastType(SRFIBL);
checkHisto(SRFIBL);
printf(" SRFIGL");
GetInputRaster(SRFIGL,nlins,ncols,stype$);
checkHisto(SRFIGL);
printf(" SRFIRL");
GetInputRaster(SRFIRL,nlins,ncols,stype$);
checkHisto(SRFIRL);
printf(" SRFINA");
GetInputRaster(SRFINA,nlins,ncols,stype$);
checkHisto(SRFINA);
printf("\n\n");
# ------------------------------------------------------------
# IF m == 2, OPEN WATER-MASK RASTER:
if (m == 2) then begin
printf("OPEN WATERMASK RASTER.");
GetInputRaster(WATERMASK,nlins,ncols,"binary");
checkHisto(WATERMASK);
IgnoreNull(WATERMASK);
printf("\n\n");
end
# ------------------------------------------------------------
# IF m == 1, CREATE NEW OUTPUT RASTER(S:
printf("CREATE NEW OUTPUT RASTER(S):\n");
if (m == 1) then begin
printf(" WATERMASK");
GetOutputRaster(WATERMASK,nlins,ncols,"binary");
IgnoreNull(WATERMASK);
CopySubobjects(SRFIBL,WATERMASK,"GEOREF");
printf(" IMAGE");
GetOutputRaster(IMAGE,nlins,ncols,"24-bit color RGB");
IgnoreNull(IMAGE);
CopySubobjects(SRFIBL,IMAGE,"GEOREF");
end
# ------------------------------------------------------------
# IF m == 2, CREATE NEW IMAGE RASTER ONLY:
if (m == 2) then begin
printf(" IMAGE");
GetOutputRaster(IMAGE,nlins,ncols,"24-bit color RGB");
IgnoreNull(IMAGE);
CopySubobjects(SRFIBL,IMAGE,"GEOREF");
end
printf("\n\n");
# ------------------------------------------------------------
# PRODUCE VALUES FOR OUTPUT WATERMASK PIXELS:
if (m == 1) then begin
printf("FILLING WATERMASK RASTER WITH VALUES.");
for each SRFINA begin
x = SRFIBL; nc = IsNull(x);
if (nc == 0) then begin
wm = 1; yt = round(facBL * x + offNA);
y = SRFINA;
if (y > yt) then wm = 0;
end
WATERMASK = wm;
end
printf("\n\n");
end
# ------------------------------------------------------------
# PRODUCE VALUES FOR OUTPUT IMAGE RASTER WATER-PIXELS:
srfiGLmaxW = 1500 / brighterW;
srfiRLmaxW = srfiGLmaxW / (redderW^1.51);
srfiBLmaxW = srfiGLmaxW * redderW;
pW = 1 / brightenDarkW;
printf("DERIVED ENHANCEMENT FACTORS FOR WATER PIXELS:\n");
printf(" srfiBLmaxW = %4d\n",srfiBLmaxW);
printf(" srfiGLmaxW = %4d\n",srfiGLmaxW);
printf(" srfiRLmaxW = %4d\n",srfiRLmaxW);
printf(" pW = %4.2f\n\n",pW);
for i=1 to 65535 begin
v = 1 + round(254*(((i - 1) / srfiBLmaxW)^pW));
if (v < 1) then v = 1;
if (v > 255) then v = 255;
sBL_8[i] = v;
end
for i=1 to 65535 begin
v = 1 + round(254*(((i - 1) / srfiGLmaxW)^pW));
if (v < 1) then v = 1;
if (v > 255) then v = 255;
sGL_8[i] = v;
end
for i=1 to 65535 begin
v = 1 + round(254*(((i - 1) / srfiRLmaxW)^pW));
if (v < 1) then v = 1;
if (v > 255) then v = 255;
sRL_8[i] = v;
end
if (tLP == 1) then begin
srfiRLmaxL = 5500 / brighterL;
srfiNAmaxL = srfiRLmaxL / redderL;
srfiGLmaxL = srfiRLmaxL * (redderL^1.01);
srfiBLmaxL = srfiRLmaxL * (redderL^1.68);
end
if (tLP == 2) then begin
srfiGLmaxL = 3951 / brighterL;
srfiBLmaxL = srfiGLmaxL * redderL;
srfiRLmaxL = srfiGLmaxL / (redderL^1.65);
srfiNAmaxL = srfiGLmaxL / (redderL^1.97);
end
if (tLP == 3) then begin
srfiGLmaxL = 3951 / brighterL;
srfiBLmaxL = srfiGLmaxL * redderL;
srfiRLmaxL = srfiGLmaxL / (redderL^1.65);
srfiNAmaxL = srfiGLmaxL / (redderL^1.97);
end
if (tLP == 4) then begin
srfiNAmaxL = 5800 / brighterL;
srfiRLmaxL = srfiNAmaxL * redderL;
srfiGLmaxL = srfiNAmaxL * (redderL^1.99);
srfiBLmaxL = srfiNAmaxL * (redderL^2.64);
end
pL = 1 / brightenDarkL;
printf("DERIVED ENHANCEMENT FACTORS FOR LAND PIXELS:\n");
if (tLP == 2) then begin
printf(" srfiBLmaxL = %4d\n",srfiBLmaxL);
end
if (tLP < 4) then begin
printf(" srfiGLmaxL = %4d\n",srfiGLmaxL);
end
if (tLP < 3) then begin
printf(" srfiRLmaxL = %4d\n",srfiRLmaxL);
end
if (tLP == 1 or tLP == 4) then begin
printf(" srfiNAmaxL = %4d\n",srfiNAmaxL);
end
printf(" pL = %4.2f\n\n",pL);
printf("COLORIZING THE WATER PIXELS.");
for each WATERMASK begin
wm = WATERMASK;
if (wm == 1) then begin
v = SRFIBL; iB = sBL_8[v]; v = SRFIGL;
iG = sGL_8[v]; v = SRFIRL; iR= sRL_8[v];
IMAGE = iR * fR + iG * fG + iB * fB;
end
end
printf("\n\n");
# ------------------------------------------------------------
# PRODUCE VALUES FOR OUTPUT IMAGE RASTERS LAND-PIXELS:
for i=1 to 65535 begin
v = 1 + round(254*(((i - 1) / srfiBLmaxL)^pL));
if (v < 1) then v = 1;
if (v > 255) then v = 255;
sBL_8[i] = v;
end
for i=1 to 65535 begin
v = 20 + round(235*(((i - 1) / srfiGLmaxL)^pL));
if (v < 1) then v = 1;
if (v > 255) then v = 255;
sGL_8[i] = v;
end
for i=1 to 65535 begin
v = 1 + round(254*(((i - 1) / srfiRLmaxL)^pL));
if (v < 1) then v = 1;
if (v > 255) then v = 255;
sRL_8[i] = v;
end
for i=1 to 65535 begin
v = 1 + round(254*(((i - 1) / srfiNAmaxL)^pL));
if (v < 1) then v = 1;
if (v > 255) then v = 255;
sNA_8[i] = v;
end
printf("COLORIZING THE LAND PIXELS.");
for each WATERMASK begin
wm = WATERMASK;
if (wm == 0) then begin
v = SRFIGL; nc = IsNull(v);
if (nc == 0) then begin
if (tLP <> 2) then begin
iB = sGL_8[v]; v = SRFIRL; iG = sRL_8[v];
v = SRFINA; iR = sNA_8[v];
end
if (tLP == 2) then begin
v = SRFIBL; iB = sBL_8[v]; v = SRFIGL;
iG = sGL_8[v]; v = SRFIRL; iR = sRL_8[v];
end
if (tLP == 3) then begin
iG = iB; iR = iB;
end
if (tLP == 4) then begin
iB = iR; iG = iR;
end
IMAGE = iR * fR + iG * fG + iB * fB;
end
end
end
printf("\n\n");
printf("CREATING HISTOGRAMS & PYRAMIDS FOR IMAGE RASTERS.");
# ------------------------------------------------------------
# CREATE HISTOGRAMS & PYRAMIDS FOR OUTPUT RASTERS:
SetNull(IMAGE,0); CreateHistogram(IMAGE,0);
CreatePyramid(IMAGE,0);
if (m == 1) then begin
CreateHistogram(WATERMASK,0); CreatePyramid(WATERMASK,0);
end
# ------------------------------------------------------------
# CLOSE or DELETE RASTERS:
CloseRaster(SRFIBL); CloseRaster(SRFIGL);
CloseRaster(SRFIRL); CloseRaster(SRFINA);
CloseRaster(WATERMASK); CloseRaster(IMAGE);
printf("\n\n");
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