PROGRAM AVRECGEN
c ***************************************************************
c * Generate scanline and pixel lookup images using AVIRIS *
c * engineering data for input into AVRECTFY (to rectify *
c * AVIRIS single band and cube data *
c * *
c * set band to band number to copy single band from BIL file *
c * band ignored when input file is a single band,fcb(13)=1 *
c * *
c * To change flightline array length, change lines: *
c * INTEGER*4 slpos2(1014,6544) and pxpos(1014,6554) *
c * maxsl = 6144 + (2*offset) *
c * This size is set for 12 (512sl) AVIRIS Segments for now *
c * *
c * K. Eric Livo Sept. 1997 *
c ***************************************************************
INTEGER*4 pixoff(256), sloff(256)
INTEGER*2 pixbuf(8192), slbuf(8192)
INTEGER*2 slpost(8192), pxpost(8192)
INTEGER*4 slposfcb(256), pxposfcb(256)
INTEGER*2 slpos(1014,6744), slpos2(1014,6744), pxpos(1014,6744)
INTEGER*4 maxsl, maxpix, cursl, curpix, offset
INTEGER*4 isl, ipix
INTEGER*4 topsl, botsl
INTEGER*4 leftpx, rghtpx, tmppx
INTEGER*4 curpxout, curslout
INTEGER*4 impixmin, impixmax, imslmin, imslmax
INTEGER*4 ifound, xdist, ydist
REAL*4 rsl, rpix
DATA pixoff,sloff/3*0,2*1,251*0,3*0,2*1,251*0/
c oversize top, left, & right sides 200 pixels
c oversize bottom 400 scanlines (output image grows downward)
offset = 200
maxsl = 6144 + (3*offset)
maxpix = 614 + (2*offset)
impixmin = 0
impixmax = 0
imslmin = 0
imslmax = 0
CALL REMLOG('=========== AVRECGEN ===========')
c sloff(1)=0 open read only (pc: to implement);
c sloff(1)=-1 open read/write
c assumes scanline and pixel offset images enter a single band each
c open scanline offset image
sloff(1) = 0
write(*,*) 'Enter SCANLINE OFFSET FILENAME'
CALL DISKIO(0,11,slbuf,sloff)
IF(sloff(1).LT.0) CALL EXIT(3)
c open pixel offset image
pixoff(1) = 0
write(*,*) 'Enter PIXEL OFFSET FILENAME'
CALL DISKIO(0,12,pixbuf,pixoff)
IF(pixoff(1).LT.0) CALL EXIT(3)
C *****
c assuming scanline and pixel-offset images are the same dimensions
c for now, just register single band inputs, not testing for BIL input
C *****
c read inputs
cursl = 0
1 CALL DISKIO(9,11,slbuf,sloff)
IF(sloff(1).EQ.-1) GOTO 2
IF(sloff(1).LT.-1) CALL EXIT(3)
CALL UNPACK(slbuf,sloff)
CALL DISKIO(9,12,pixbuf,pixoff)
IF(pixoff(1).EQ.-1) GOTO 2
IF(pixoff(1).LT.-1) CALL EXIT(3)
CALL UNPACK(pixbuf,pixoff)
cursl = cursl + 1
c process inputs
c X & Y projection offsets are truncated to 'OFFSET'(200 pixel) maximum
do 10,curpix=1,sloff(2)
rsl = (slbuf(curpix) + 0.5) / 10.0
isl = int(rsl)
curslout = cursl + offset + isl
if(curslout .lt. 1) then
write(*,*) 'output scanline projected above top of image',
& curslout
curslout = 1
elseif(curslout .gt. maxsl) then
write(*,*) 'output scanline projected below bottom of',
& ' image', curslout
curslout = maxsl
endif
rpix = (pixbuf(curpix) + 0.5) / 10.0
ipix = int(rpix)
curpxout = curpix+offset+ipix
if(curpxout .lt. 1 ) then
write(*,*) 'output pixel projected to left of image',
& curpxout
curpxout = 1
elseif(curpxout .gt. maxpix) then
write(*,*) 'output pixel projected to right of image',
& curpxout
curpxout = maxpix
endif
c set 1st data point min and max
if (impixmin.eq.0) impixmin=curpxout
if (impixmax.eq.0) impixmax=curpxout
if (imslmin.eq.0) imslmin=curslout
if (imslmax.eq.0) imslmax=curslout
c find min and max pix's and sl's
if (curpxout.lt.impixmin) impixmin=curpxout
if (curpxout.gt.impixmax) impixmax=curpxout
if (curslout.lt.imslmin) imslmin=curslout
if (curslout.gt.imslmax) imslmax=curslout
c assign original-sl, and original-sample
slpos(curpxout, curslout) = cursl
slpos2(curpxout, curslout) = cursl
pxpos(curpxout, curslout) = curpix
10 continue
goto 1
2 continue
c open output: scanline lookup (sl's of orig unrectified image)
slposfcb(1)=16
slposfcb(2)=(impixmax-impixmin+1)
slposfcb(13)=1
write(*,*) 'Enter scanline LOOKUP output image'
CALL DISKIO(0,14,slpost,slposfcb)
IF(slposfcb(1).LT.0) CALL EXIT(3)
c open output: pixel lookup (pix's of orig unrectified image)
pxposfcb(1)=16
pxposfcb(2)=(impixmax-impixmin+1)
pxposfcb(13)=1
write(*,*) 'Enter pixel LOOKUP output image'
CALL DISKIO(0,15,pxpost,pxposfcb)
IF(pxposfcb(1).LT.0) CALL EXIT(3)
c find DN 0's line top and bottom (1st non-zero value)
c and left and right
c note: slpos2(pix,sl) = 0 only where holes are formed through warping
c original image to rectified output image
do 20, cursl=imslmin,imslmax
do 30, curpix=impixmin,impixmax
if(slpos2(curpix,cursl) .eq. 0) then
topsl = cursl
ifound = 0
do while((ifound.eq.0).and.(topsl.gt.imslmin))
topsl = topsl - 1
ifound = slpos2(curpix,topsl)
enddo
if(ifound.eq.0) topsl = 0
botsl = cursl
ifound = 0
do while((ifound.eq.0).and.(botsl.lt.imslmax))
botsl = botsl + 1
ifound = slpos2(curpix,botsl)
enddo
if(ifound.eq.0) botsl = 0
leftpx = curpix
ifound = 0
do while((ifound.eq.0).and.(leftpx.gt.impixmin))
leftpx = leftpx - 1
ifound = slpos2(leftpx,cursl)
enddo
if(ifound.eq.0) leftpx = 0
rghtpx = curpix
ifound = 0
do while((ifound.eq.0).and.(rghtpx.lt.impixmax))
rghtpx = rghtpx + 1
ifound = slpos2(rghtpx,cursl)
enddo
if(ifound.eq.0) rghtpx = 0
ydist = 0
if((topsl.gt.0).and.(botsl.gt.0)) then
ydist = topsl - cursl
if(IABS(ydist).gt.(botsl-cursl)) ydist=botsl-cursl
endif
xdist = 0
if((leftpx.gt.0).and.(rghtpx.gt.0)) then
xdist = leftpx - curpix
if(IABS(xdist).gt.(rghtpx-curpix)) xdist=rghtpx-curpix
endif
if((xdist.ne.0).and.(ydist.ne.0)) then
if(IABS(xdist) .le. IABS(ydist)) then
slpos(curpix,cursl) = slpos2(curpix+xdist,cursl)
pxpos(curpix,cursl) = pxpos(curpix+xdist,cursl)
else
slpos(curpix,cursl) = slpos2(curpix,cursl+ydist)
pxpos(curpix,cursl) = pxpos(curpix,cursl+ydist)
endif
else
if(xdist.ne.0) then
slpos(curpix,cursl) = slpos2(curpix+xdist,cursl)
pxpos(curpix,cursl) = pxpos(curpix+xdist,cursl)
endif
if(ydist.ne.0) then
slpos(curpix,cursl) = slpos2(curpix,cursl+ydist)
pxpos(curpix,cursl) = pxpos(curpix,cursl+ydist)
endif
endif
c NOTE: if no match, slpos remains set to sl=0 & pxpos set to pix=0
endif
30 continue
20 continue
c write registered output
do 200, cursl=imslmin,imslmax
do 210, curpix=impixmin,impixmax
tmppx = curpix - impixmin + 1
slpost(tmppx) = slpos(curpix,cursl)
pxpost(tmppx) = pxpos(curpix,cursl)
210 continue
CALL DISKIO(10,14,slpost,slposfcb)
IF(slposfcb(1).LT.0) CALL EXIT(4)
CALL DISKIO(10,15,pxpost,pxposfcb)
IF(pxposfcb(1).LT.0) CALL EXIT(4)
200 continue
c close all files
CALL DISKIO(6,11,slbuf,sloff)
CALL DISKIO(6,12,pixbuf,pixoff)
CALL DISKIO(6,14,slpost,slposfcb)
CALL DISKIO(6,15,pxpost,pxposfcb)
STOP
END
U.S. Geological Survey,
a bureau of the U.S. Department of the Interior
This page URL= http://speclab.cr.usgs.gov/aviris.geom.software/avrecgen.f.html
This page is maintained by: Eric Livo elivo@speclab.cr.usgs.gov
Last modified December 16, 1998.