include	<ctype.h>
include	<foc.h>

# RINTERX -- performs linear interpolation between reseau grids
#
# D. Giaretta, 01-Aug-1987	Original SPP version

procedure t_rinterx()

char	input1[SZ_FNAME]		# input file 1
char	entry1[SZ_R_PATTERN]		# input1 entry template
char	input2[SZ_FNAME]		# input file 
char	entry2[SZ_R_PATTERN]		# input2 entry template
char	output[SZ_FNAME]		# output file 
char	outentry[SZ_R_PATTERN]		# input2 entry template
real	mag1, mag2			# mag. factors along x,y
real	shift1, shift2			# shifts in x and y
#--

pointer	rpout, rp1, rp2			# res. structure pointers
pointer rentp1, rentp2			# res. entry structure
int	nrows, ncols			# rows, cols in grid
bool	streq()
int	nowhite(), i
real	clgetr()
bool	rs_check()
pointer	rs_open(), rs_template(), rlist
int	rs_readr(), rs_list()

begin

	call clgstr( "input1", input1, SZ_FNAME)
	call clgstr( "entry1", entry1, SZ_R_PATTERN)
	call clgstr( "input2", input2, SZ_FNAME)
	call clgstr( "entry2", entry2, SZ_R_PATTERN)

	call clgstr( "output", output, SZ_FNAME)
	call clgstr( "outentry", outentry, SZ_R_PATTERN)
	# default output entry name is the input1 entry name
	if ( nowhite( outentry, outentry, SZ_R_PATTERN) == 0 )
	    call strcpy( entry1, outentry, SZ_R_PATTERN)

	mag1   = clgetr("mag1")
	mag2   = clgetr("mag2")
	shift1 = clgetr("shift1")
	shift2 = clgetr("shift2")
	
	if ( streq( input1, output) ) {
	    rp1 = rs_open( input1, READ_WRITE, 0)
	    if ( streq( input1, input2) )
		rp2 = rp1
	    rpout = rp1
	} else if ( streq( input2, output) ) {
	    rp2 = rs_open( input2, READ_WRITE, 0)
	    rpout = rp2
	} else {
	    rp1   = rs_open( input1,  READ_ONLY, 0)
	    if ( streq( input1, input2) )
		rp2 = rp1
	    else
		rp2   = rs_open( input2,  READ_ONLY, 0)
	    rpout = rs_open( output, NEW_COPY, rp1)
	}
	
	# check grid sizes match
	nrows = RES_NROWS(rp1)
	ncols = RES_NCOLS(rp1)
	if ( !rs_check( rp2, RES_NROWS(rp1), RES_NCOLS(rp1) ) )
	    call error( 0, "reseau grids must have same size")

	# get space for input and output coords
	call rs_alloc( rp1, TY_REAL, rentp1)
	call rs_alloc( rp2, TY_REAL, rentp2)

	rlist = rs_template( rp2, entry2)
	if ( rs_readr( rp2, rentp2, EOS, rs_list(rlist) ) == RES_F_NONEXTENTRY )
	    call error ( 0, "cannot read entry2 of infile2")

	rlist = rs_template( rp1, entry1)
	if ( rs_readr( rp1, rentp1, EOS, rs_list(rlist) ) == RES_F_NONEXTENTRY )
	    call error ( 0, "cannot read entry1 of infile1")

	# perform the calculations
	do i = 1, nrows*ncols {
	    if ( !Memb[RES_FPT(rentp1)+i-1] && !Memb[RES_FPT(rentp2)+i-1] ) {
		Memr[RES_XPT(rentp1)+i-1] = mag1*Memr[RES_XPT(rentp1)+i-1] +
					mag2*Memr[RES_XPT(rentp2)+i-1] + shift1
		Memr[RES_YPT(rentp1)+i-1] = mag1*Memr[RES_YPT(rentp1)+i-1] +
					mag2*Memr[RES_YPT(rentp2)+i-1] + shift2
	    } else {
		Memb[RES_FPT(rentp1)+i-1] = true
		Memr[RES_XPT(rentp1)+i-1] = INDEFR
		Memr[RES_YPT(rentp1)+i-1] = INDEFR
	    }
	}

	# update the tracking info
	call rs_t_update(rentp1, "I", SZ_R_TRACKING)
	# insert name
	call strcpy( outentry, RES_ENT_ENTRY(rentp1), SZ_R_ENTRY)
	# write entry out
	call rs_write ( rpout, rentp1, 0)

	# close the files
	call rs_close(rpout)
	call rs_close(rp1)
	call rs_close(rp2)

end