aboutsummaryrefslogtreecommitdiff
path: root/pkg/utilities/nttools/tmatch
diff options
context:
space:
mode:
Diffstat (limited to 'pkg/utilities/nttools/tmatch')
-rw-r--r--pkg/utilities/nttools/tmatch/getmatch.x101
-rw-r--r--pkg/utilities/nttools/tmatch/getnorm.x67
-rw-r--r--pkg/utilities/nttools/tmatch/getweight.x96
-rw-r--r--pkg/utilities/nttools/tmatch/infomatch.x219
-rw-r--r--pkg/utilities/nttools/tmatch/mkpkg20
-rw-r--r--pkg/utilities/nttools/tmatch/putmatch.x102
-rw-r--r--pkg/utilities/nttools/tmatch/rowname.x61
-rw-r--r--pkg/utilities/nttools/tmatch/setindex.x13
-rw-r--r--pkg/utilities/nttools/tmatch/sortclose.x50
-rw-r--r--pkg/utilities/nttools/tmatch/sortdist.x50
-rw-r--r--pkg/utilities/nttools/tmatch/tmatch.x138
11 files changed, 917 insertions, 0 deletions
diff --git a/pkg/utilities/nttools/tmatch/getmatch.x b/pkg/utilities/nttools/tmatch/getmatch.x
new file mode 100644
index 00000000..083d4156
--- /dev/null
+++ b/pkg/utilities/nttools/tmatch/getmatch.x
@@ -0,0 +1,101 @@
+#* HISTORY *
+#* B.Simon 24-Aug-94 Original
+#* B.Simon 18-Sep-00 Revised search termination criterion
+
+# GETMATCH -- Find rows in second table witch match rows in first
+
+procedure getmatch (in1, in2, ncol, col1, col2, weight, nrow1, index1,
+ nrow2, index2, maxnorm, sphere, closest, dist)
+
+pointer in1 # i: first table descriptor
+pointer in2 # i: second table descriptor
+int ncol # i: number of match columns
+pointer col1[ARB] # i: match columns in first table
+pointer col2[ARB] # i: match columns in second table
+double weight[ARB] # i: weights used in computing norm
+int nrow1 # i: number of rows in first table
+int index1[ARB] # i: sorted row indices for first table
+int nrow2 # i: number of rows in second table
+int index2[ARB] # i: sorted row indices for second table
+double maxnorm # i: maximum norm used in match
+bool sphere # i: apply spherical correction to first column?
+int closest[ARB] # o: closest match in second table to first
+double dist[ARB] # o: distance between matched rows
+#--
+double sqnorm, proj, abnorm, norm, minabnorm, minnorm
+int idx, jdx, irow, jrow, krow, jlast
+
+begin
+ jlast = 1
+ sqnorm = maxnorm * maxnorm
+
+ # Find the row in the second table which minimizes the norm for
+ # each row of the first table
+
+
+ do idx = 1, nrow1 {
+ irow = index1[idx]
+ jrow = index2[jlast]
+
+ # The initial guess is the row which matched last time
+
+ call getnorm (in1, in2, ncol, col1, col2, irow, jrow,
+ weight, sphere, proj, abnorm, norm)
+
+ minabnorm = abnorm
+ minnorm = norm
+ krow = jrow
+
+ # Search backwards for a row which minimizes the norm
+ # Terminate the search when the first dimension of the norm (proj)
+ # is greater than the minimum norm, as all subsequent rows have
+ # norms that must be greater than the minimum
+
+ do jdx = jlast-1, 1, -1 {
+ jrow = index2[jdx]
+
+ call getnorm (in1, in2, ncol, col1, col2, irow, jrow,
+ weight, sphere, proj, abnorm, norm)
+
+ if (proj > minabnorm)
+ break
+
+ if (norm < minnorm) {
+ minabnorm = abnorm
+ minnorm = norm
+ krow = jrow
+ }
+ }
+
+ # Search forwards for a row that minimizes the norm
+ # Use the same termination condition as for the forwards search
+
+ do jdx = jlast+1, nrow2 {
+ jrow = index2[jdx]
+
+ call getnorm (in1, in2, ncol, col1, col2, irow, jrow,
+ weight, sphere, proj, abnorm, norm)
+
+ if (proj > minabnorm)
+ break
+
+ if (norm < minnorm) {
+ minabnorm = abnorm
+ minnorm = norm
+ krow = jrow
+ }
+ }
+
+ if (minnorm > sqnorm) {
+ dist[irow] = maxnorm
+ closest[irow] = 0
+ } else {
+ dist[irow] = sqrt (minnorm)
+ closest[irow] = krow
+ }
+
+ jlast = krow
+ }
+
+end
+
diff --git a/pkg/utilities/nttools/tmatch/getnorm.x b/pkg/utilities/nttools/tmatch/getnorm.x
new file mode 100644
index 00000000..fed277c6
--- /dev/null
+++ b/pkg/utilities/nttools/tmatch/getnorm.x
@@ -0,0 +1,67 @@
+include <math.h>
+
+#* HISTORY *
+#* B.Simon 24-Aug-94 original
+#* B.Simon 18-Sep-00 Revised computation of proj and added abnorm
+
+# GETNORM -- Compute the squared norm between two table rows
+
+procedure getnorm (in1, in2, ncol, col1, col2, row1, row2, weight, sphere,
+ proj, abnorm, norm)
+
+pointer in1 # i: first table descriptor
+pointer in2 # i: second table descriptor
+int ncol # i: number of match columns
+pointer col1[ARB] # i: match columns in first table
+pointer col2[ARB] # i: match columns in second table
+int row1 # i: row number in first table
+int row2 # i: row number in second table
+double weight[ARB] # i: weights used in computing norm
+bool sphere # i: apply spherical correction to first column?
+double proj # o: projection of norm on first axis
+double abnorm # o: norm, possibly without spherical correction
+double norm # o: norm (distance) between rows in two tables
+#--
+int i
+double val1, val2, dif
+
+begin
+ # Calculate first component of norm
+
+ call tbegtd (in1, col1[1], row1, val1)
+ call tbegtd (in2, col2[1], row2, val2)
+
+ dif = weight[1] * (val1 - val2)
+ proj = dif * dif
+ abnorm = proj
+
+ # Apply correction for spherical coordinates
+
+ if (sphere) {
+ if (dif > 180) {
+ dif = dif - 360
+ } else if (dif < -180) {
+ dif = dif + 360
+ }
+
+ call tbegtd (in1, col1[2], row1, val1)
+ call tbegtd (in2, col2[2], row2, val2)
+
+ val1 = 0.5 * weight[2] * (val1 + val2)
+ dif = dif * cos (DEGTORAD(val1))
+ }
+
+ # Compute remaining components
+
+ norm = dif * dif
+ do i = 2, ncol {
+ call tbegtd (in1, col1[i], row1, val1)
+ call tbegtd (in2, col2[i], row2, val2)
+
+ dif = weight[i] * (val1 - val2)
+ abnorm = abnorm + dif * dif
+ norm = norm + dif * dif
+ }
+
+end
+
diff --git a/pkg/utilities/nttools/tmatch/getweight.x b/pkg/utilities/nttools/tmatch/getweight.x
new file mode 100644
index 00000000..ef8b9204
--- /dev/null
+++ b/pkg/utilities/nttools/tmatch/getweight.x
@@ -0,0 +1,96 @@
+include <math.h>
+include <tbset.h>
+
+#* HISTORY *
+#* B.Simon 24-Aug-94 original
+
+# GETWEIGHT -- Get array of weights from list of factors or tables
+
+procedure getweight (ncol, col1, col2, factor, weight)
+
+int ncol # i: number of match columns
+pointer col1[ARB] # i: match columns from first table
+pointer col2[ARB] # i: match columns from second table
+char factor[ARB] # i: list of factors
+double weight[ARB] # o: array of weights
+#--
+double unitval[6]
+int invert[6]
+int ic, jc, nc, icol, jcol, type1, type2, item
+pointer sp, value, unit1, unit2, errmsg
+
+
+data unitval / 1.0, 3600.0, 60.0, 1.0, 15.0, RADIAN /
+data invert / NO, YES, YES, NO, NO, NO /
+
+string unitlist "|seconds|minutes|degrees|hours|radians|"
+string badvalue "Value in factor string is not a number (%s)"
+string badunits "Units mismatch in column %d of tables"
+
+int ctod(), word_fetch(), strdic()
+
+begin
+ # Allocate memory for temporary strings
+
+ call smark (sp)
+ call salloc (value, SZ_FNAME, TY_CHAR)
+ call salloc (unit1, SZ_FNAME, TY_CHAR)
+ call salloc (unit2, SZ_FNAME, TY_CHAR)
+ call salloc (errmsg, SZ_LINE, TY_CHAR)
+
+ # Get each string from the list and convert to a number
+
+ ic = 1
+ icol = 0
+ while (word_fetch (factor, ic, Memc[value], SZ_FNAME) > 0) {
+ icol = icol + 1
+
+ jc = 1
+ nc = ctod (Memc[value], jc, weight[icol])
+ if (Memc[value+jc-1] != EOS) {
+ call sprintf (Memc[errmsg], SZ_LINE, badvalue)
+ call pargstr (Memc[value])
+
+ call error (1, Memc[errmsg])
+ }
+ }
+
+ # Set remaining weights according to column units
+
+ do jcol = icol+1, ncol {
+ # Read units from table
+
+ call tbcigt (col1[jcol], TBL_COL_UNITS, Memc[unit1], SZ_FNAME)
+ call tbcigt (col2[jcol], TBL_COL_UNITS, Memc[unit2], SZ_FNAME)
+
+ # Search for units in dictionary
+
+ call strlwr (Memc[unit1])
+ call strlwr (Memc[unit2])
+
+ type1 = strdic (Memc[unit1], Memc[unit1], SZ_FNAME, unitlist)
+ type2 = strdic (Memc[unit2], Memc[unit2], SZ_FNAME, unitlist)
+
+ # Take exit if units do not match
+
+ if (type1 != type2) {
+ call sprintf (Memc[errmsg], SZ_LINE, badunits)
+ call pargi (jcol)
+
+ call error (1, Memc[errmsg])
+ }
+
+ # Read corresponding weight from unit value array
+ # The first weight (1.0) is for missing or unknown units
+
+ item = type1 + 1
+
+ if (invert[item] == NO) {
+ weight[jcol] = unitval[item]
+ } else {
+ weight[jcol] = 1.0 / unitval[item]
+ }
+ }
+
+ call sfree (sp)
+end
diff --git a/pkg/utilities/nttools/tmatch/infomatch.x b/pkg/utilities/nttools/tmatch/infomatch.x
new file mode 100644
index 00000000..d90f89a8
--- /dev/null
+++ b/pkg/utilities/nttools/tmatch/infomatch.x
@@ -0,0 +1,219 @@
+include <tbset.h>
+
+#* HISTORY *
+#* B.Simon 25-Aug-94 original
+
+# INFOMATCH -- Print diagnostic information for tmatch
+
+procedure infomatch (diagfile, in1, in2, nmcol1, nmcol2, maxnorm,
+ nclosest, closest, dist)
+
+char diagfile[ARB] # i: diagnostic output file
+pointer in1 # i: first table's descriptor
+pointer in2 # i: second table's descriptor
+char nmcol1[ARB] # i: name columns in first table
+char nmcol2[ARB] # i: name columns in second table
+double maxnorm # i: maximum allowed distance between matched rows
+int nclosest # i: length of closest array
+int closest[ARB] # i: array of closest matches between tables
+double dist[ARB] # i: distance between matched rows
+#--
+bool first, same
+int fd, namelen, mxcol1, mxcol2, ncol1, ncol2, idx, jdx, irow, jrow
+pointer sp, index, name, col1, col2
+
+string ziptitle "\nThe following objects were not matched:\n"
+string duptitle "\nThe following objects matched the same object:\n"
+string bigtitle "\nThe following objects have the largest norms:\n"
+string normfmt "Norm = %0.7g\n"
+string rowformat "%d:%d %s\n"
+
+bool is_blank()
+int open(), envgeti(), tbpsta()
+
+begin
+ # Open the diagnostics file
+
+ if (is_blank (diagfile))
+ return
+
+ fd = open (diagfile, WRITE_ONLY, TEXT_FILE)
+
+ # Get maximum length of diagnostic string
+
+ iferr {
+ namelen = envgeti ("ttyncols") - 10
+ } then {
+ namelen = 70
+ }
+
+ # Allocate dynamic memory
+
+ call smark (sp)
+ call salloc (index, nclosest, TY_INT)
+ call salloc (name, namelen, TY_CHAR)
+
+ # Get column descriptors for name columns
+
+ mxcol1 = tbpsta (in1, TBL_NCOLS)
+ mxcol2 = tbpsta (in2, TBL_NCOLS)
+
+ call salloc (col1, mxcol1, TY_INT)
+ call salloc (col2, mxcol2, TY_INT)
+
+ if (is_blank (nmcol1)) {
+ ncol1 = 0
+ } else {
+ call tctexp (in1, nmcol1, mxcol1, ncol1, Memi[col1])
+ }
+
+ if (is_blank (nmcol2)) {
+ ncol2 = 0
+ } else {
+ call tctexp (in2, nmcol2, mxcol2, ncol2, Memi[col2])
+ }
+
+ # Sort the closest array
+
+ call setindex (Memi[index], nclosest)
+ call sortclose (nclosest, closest, Memi[index])
+
+ # Print the objects that were not matched
+
+ first = true
+ do idx = 1, nclosest {
+ irow = Memi[index+idx-1]
+ if (closest[irow] != 0)
+ break
+
+ if (first) {
+ first = false
+ call fprintf (fd, ziptitle)
+ }
+
+ call rowname (in1, irow, ncol1, Memi[col1], Memc[name], namelen)
+ call fprintf (fd, rowformat)
+ call pargi (1)
+ call pargi (irow)
+ call pargstr (Memc[name])
+ }
+
+ # Print the objects which are matched more than once
+
+ same = false
+ first = true
+ do idx = 2, nclosest {
+ irow = Memi[index+idx-1]
+ jrow = Memi[index+idx-2]
+
+ if (closest[irow] == 0)
+ next
+
+ if (closest[irow] == closest[jrow]) {
+ same = true
+
+ if (first) {
+ first = false
+ call fprintf (fd, duptitle)
+ }
+
+ call rowname (in1, jrow, ncol1, Memi[col1],
+ Memc[name], namelen)
+
+ call fprintf (fd, rowformat)
+ call pargi (1)
+ call pargi (jrow)
+ call pargstr (Memc[name])
+
+ } else if (same) {
+ same = false
+
+ call rowname (in1, jrow, ncol1, Memi[col1],
+ Memc[name], namelen)
+
+ call fprintf (fd, rowformat)
+ call pargi (1)
+ call pargi (jrow)
+ call pargstr (Memc[name])
+
+ call rowname (in2, closest[jrow], ncol2, Memi[col2],
+ Memc[name], namelen)
+
+ call fprintf (fd, rowformat)
+ call pargi (2)
+ call pargi (closest[jrow])
+ call pargstr (Memc[name])
+
+ call fprintf (fd, "\n")
+ }
+ }
+
+ if (same) {
+ same = false
+ irow = Memi[index+nclosest-1]
+
+ call rowname (in1, irow, ncol1, Memi[col1],
+ Memc[name], namelen)
+
+ call fprintf (fd, rowformat)
+ call pargi (1)
+ call pargi (irow)
+ call pargstr (Memc[name])
+
+ call rowname (in2, closest[irow], ncol2, Memi[col2],
+ Memc[name], namelen)
+
+ call fprintf (fd, rowformat)
+ call pargi (2)
+ call pargi (closest[irow])
+ call pargstr (Memc[name])
+
+ call fprintf (fd, "\n")
+ }
+
+ # Sort the dist array
+
+ call setindex (Memi[index], nclosest)
+ call sortdist (nclosest, dist, Memi[index])
+
+ # Print the ten objects with the largest norms
+
+ jdx = 0
+ do idx = nclosest, 1, -1 {
+ irow = Memi[index+idx-1]
+ if (dist[irow] == maxnorm)
+ next
+
+ if (jdx == 0)
+ call fprintf (fd, bigtitle)
+
+ jdx = jdx + 1
+ if (jdx > 10)
+ break
+
+ call fprintf (fd, normfmt)
+ call pargd (dist[irow])
+
+ call rowname (in1, irow, ncol1, Memi[col1],
+ Memc[name], namelen)
+
+ call fprintf (fd, rowformat)
+ call pargi (1)
+ call pargi (irow)
+ call pargstr (Memc[name])
+
+ call rowname (in2, closest[irow], ncol2, Memi[col2],
+ Memc[name], namelen)
+
+ call fprintf (fd, rowformat)
+ call pargi (2)
+ call pargi (closest[irow])
+ call pargstr (Memc[name])
+
+ call fprintf (fd, "\n")
+
+ }
+
+ call close (fd)
+ call sfree (sp)
+end
diff --git a/pkg/utilities/nttools/tmatch/mkpkg b/pkg/utilities/nttools/tmatch/mkpkg
new file mode 100644
index 00000000..b6e6be18
--- /dev/null
+++ b/pkg/utilities/nttools/tmatch/mkpkg
@@ -0,0 +1,20 @@
+# Update the tmatch application code in the ttools package library
+# Author: B.Simon, 30-Aug-94
+
+$checkout libpkg.a ../
+$update libpkg.a
+$checkin libpkg.a ../
+$exit
+
+libpkg.a:
+ getmatch.x
+ getnorm.x <math.h>
+ getweight.x <math.h> <tbset.h>
+ infomatch.x <tbset.h>
+ putmatch.x <tbset.h>
+ rowname.x
+ setindex.x
+ sortclose.x
+ sortdist.x
+ tmatch.x <tbset.h>
+ ;
diff --git a/pkg/utilities/nttools/tmatch/putmatch.x b/pkg/utilities/nttools/tmatch/putmatch.x
new file mode 100644
index 00000000..27407107
--- /dev/null
+++ b/pkg/utilities/nttools/tmatch/putmatch.x
@@ -0,0 +1,102 @@
+include <tbset.h>
+
+#* HISTORY *
+# B.Simon 25-Aug-94 Original
+
+# PUTMATCH -- Write matched rows in input as a single row in output table
+
+procedure putmatch (output, incol1, incol2, in1, in2, nclosest, closest)
+
+char output[ARB] # i: output table name
+char incol1[ARB] # i: list of columns to copy from first table
+char incol2[ARB] # i: list of columns to copy from second table
+pointer in1 # i: first table's descriptor
+pointer in2 # i: second table's descriptor
+int nclosest # i: length of closest array
+int closest[ARB] # i: indices of rows in second table closest to first
+#--
+int mxcol1, mxcol2, maxcol, ncol1, ncol2, ncol, type1, type2
+pointer colnum, datatype, lendata, lenfmt, icol, irow, jrow
+pointer sp, colname, colunits, colfmt, oldcol, newcol,out
+
+string nomatch "WARNING: No rows matched between tables, output \
+table is empty\n"
+
+int tbpsta()
+pointer tbtopn()
+
+begin
+ # Allocate memory for temporary strings
+
+ call smark (sp)
+ call salloc (colname, SZ_COLNAME, TY_CHAR)
+ call salloc (colunits, SZ_COLUNITS, TY_CHAR)
+ call salloc (colfmt, SZ_COLFMT, TY_CHAR)
+
+ # Get column descriptors from input tables
+
+ mxcol1 = tbpsta (in1, TBL_NCOLS)
+ mxcol2 = tbpsta (in2, TBL_NCOLS)
+ maxcol = mxcol1 + mxcol2
+
+ call salloc (oldcol, maxcol, TY_INT)
+ call salloc (newcol, maxcol, TY_INT)
+
+ call tctexp (in1, incol1, mxcol1, ncol1, Memi[oldcol])
+ call tctexp (in2, incol2, mxcol2, ncol2, Memi[oldcol+ncol1])
+ ncol = ncol1 + ncol2
+
+ # Create output table
+
+ out = tbtopn (output, NEW_FILE, NULL)
+
+ # Set type (text, row ordered, column ordered)
+
+ type1 = tbpsta (in1, TBL_WHTYPE)
+ type2 = tbpsta (in2, TBL_WHTYPE)
+ if (type1 == type2)
+ call tbpset (out, TBL_WHTYPE, type1)
+
+ # Create columns in output table
+
+ do icol = 1, ncol {
+ call tbcinf (Memi[oldcol+icol-1], colnum, Memc[colname],
+ Memc[colunits], Memc[colfmt], datatype,
+ lendata, lenfmt)
+
+ call newcolnam (ncol, Memi[oldcol], icol,
+ Memc[colname], SZ_COLNAME)
+
+ call tbcdef (out, Memi[newcol+icol-1], Memc[colname],
+ Memc[colunits], Memc[colfmt], datatype, lendata, 1)
+ }
+
+ # Copy header keywords from first input table
+
+ call tbtcre (out)
+ call tbhcal (in1, out)
+
+ # Copy rows from input table to output
+
+ jrow = 0
+ do irow = 1, nclosest {
+ if (closest[irow] == 0)
+ next
+
+ jrow = jrow + 1
+ call tbrcsc (in1, out, Memi[oldcol], Memi[newcol],
+ irow, jrow, ncol1)
+ call tbrcsc (in2, out, Memi[oldcol+ncol1], Memi[newcol+ncol1],
+ closest[irow], jrow, ncol2)
+ }
+
+ # Write warning message if no rows matched
+
+ if (jrow == 0)
+ call eprintf (nomatch)
+
+ # Clean up
+
+ call tbtclo (out)
+ call sfree (sp)
+end
diff --git a/pkg/utilities/nttools/tmatch/rowname.x b/pkg/utilities/nttools/tmatch/rowname.x
new file mode 100644
index 00000000..484fcd75
--- /dev/null
+++ b/pkg/utilities/nttools/tmatch/rowname.x
@@ -0,0 +1,61 @@
+#* HISTORY *
+#* B.Simon 26-Aug-94
+
+# ROWNAME -- Create name for table row by concatenating column values
+
+procedure rowname (in, irow, ncol, col, name, namelen)
+
+pointer in # i: table descriptor
+int irow # i: table row number
+int ncol # i: number of table columns
+pointer col[ARB] # i: table column pointers
+char name[ARB] # o: concatenated values of columns
+int namelen # i: maximum name length
+#--
+int ic, jc, icol
+pointer sp, value
+
+begin
+ # Allocate memory for column buffer
+
+ call smark (sp)
+ call salloc (value, SZ_LINE, TY_CHAR)
+
+ # Concatenate column values into name string
+
+ jc = 0
+ icol = 0
+ for (ic = 1; ic <= namelen; ic = ic + 1) {
+
+ # A value of zero is a flag to read the next coumn
+
+ if (jc == 0) {
+ icol = icol + 1
+ if (icol > ncol) {
+ if (ic > 1)
+ ic = ic - 1 # remove trailing blank
+
+ break
+ }
+
+ call tbegtt (in, col[icol], irow, Memc[value], SZ_LINE)
+ }
+
+ # Copy a single character from the buffer to the output string
+ # until the buffer is exhausted. At this point we copy a blank
+ # as a spacer and set the value of jc as a flag to read the
+ # next column.
+
+ if (Memc[value+jc] == EOS) {
+ name[ic] = ' '
+ jc = 0
+ } else {
+ name[ic] = Memc[value+jc]
+ jc = jc + 1
+ }
+ }
+
+ name[ic] = EOS
+ call sfree (sp)
+end
+
diff --git a/pkg/utilities/nttools/tmatch/setindex.x b/pkg/utilities/nttools/tmatch/setindex.x
new file mode 100644
index 00000000..9640f1cd
--- /dev/null
+++ b/pkg/utilities/nttools/tmatch/setindex.x
@@ -0,0 +1,13 @@
+# SETINDEX -- Initialize an index array
+
+procedure setindex (index, len)
+
+int index[ARB] # o: index rray
+int len # i: array length
+#--
+int i
+
+begin
+ do i = 1, len
+ index[i] = i
+end
diff --git a/pkg/utilities/nttools/tmatch/sortclose.x b/pkg/utilities/nttools/tmatch/sortclose.x
new file mode 100644
index 00000000..303ac250
--- /dev/null
+++ b/pkg/utilities/nttools/tmatch/sortclose.x
@@ -0,0 +1,50 @@
+#* HISTORY *
+#* B.Simon 25-Aug-94 original
+# Phil Hodge 12-Jul-2005 in sortclose, add 'int cmpclose()'
+
+# SORTCLOSE -- Sort the closest array
+
+procedure sortclose (nclosest, closest, index)
+
+int nclosest # i: length of closest and index arrays
+int closest[ARB] # i: indices of second table rows matching first
+int index[ARB] # u: indices of first table rows, in sort order on exit
+#--
+int cmpclose()
+extern cmpclose
+pointer sp, close2
+
+begin
+ call smark (sp)
+ call salloc (close2, nclosest, TY_INT)
+
+ call amovi (closest, Memi[close2], nclosest)
+ call gqsort (index, nclosest, cmpclose, close2)
+
+ call sfree (sp)
+end
+
+# CMPCLOSE -- Compare two elements in the close array
+
+int procedure cmpclose (closest, ielem, jelem)
+
+pointer closest # address of the closest array
+int ielem # first element
+int jelem # second element
+#--
+int order
+
+begin
+ if (Memi[closest+ielem-1] < Memi[closest+jelem-1]) {
+ order = -1
+ } else if (Memi[closest+ielem-1] > Memi[closest+jelem-1]) {
+ order = 1
+ } else if (ielem < jelem) {
+ order = -1
+ } else if (ielem < jelem) {
+ order = 1
+ } else {
+ order = 0
+ }
+ return (order)
+end
diff --git a/pkg/utilities/nttools/tmatch/sortdist.x b/pkg/utilities/nttools/tmatch/sortdist.x
new file mode 100644
index 00000000..b0745676
--- /dev/null
+++ b/pkg/utilities/nttools/tmatch/sortdist.x
@@ -0,0 +1,50 @@
+#* HISTORY *
+#* B.Simon 25-Aug-94 original
+# Phil Hodge 28-Sept-2005 in sortdist, add 'int cmpdist()'
+
+# SORTDIST -- Sort the dist array
+
+procedure sortdist (ndist, dist, index)
+
+int ndist # i: length of dist and index arrays
+double dist[ARB] # i: indices of second table rows matching first
+int index[ARB] # u: indices of first table rows, in sort order on exit
+#--
+int cmpdist()
+extern cmpdist
+pointer sp, dist2
+
+begin
+ call smark (sp)
+ call salloc (dist2, ndist, TY_DOUBLE)
+
+ call amovd (dist, Memd[dist2], ndist)
+ call gqsort (index, ndist, cmpdist, dist2)
+
+ call sfree (sp)
+end
+
+# CMPDIST -- Compare two elements in the dist array
+
+int procedure cmpdist (dist, ielem, jelem)
+
+pointer dist # address of the dist array
+int ielem # first element
+int jelem # second element
+#--
+int order
+
+begin
+ if (Memd[dist+ielem-1] < Memd[dist+jelem-1]) {
+ order = -1
+ } else if (Memd[dist+ielem-1] > Memd[dist+jelem-1]) {
+ order = 1
+ } else if (ielem < jelem) {
+ order = -1
+ } else if (ielem < jelem) {
+ order = 1
+ } else {
+ order = 0
+ }
+ return (order)
+end
diff --git a/pkg/utilities/nttools/tmatch/tmatch.x b/pkg/utilities/nttools/tmatch/tmatch.x
new file mode 100644
index 00000000..2d3ea22d
--- /dev/null
+++ b/pkg/utilities/nttools/tmatch/tmatch.x
@@ -0,0 +1,138 @@
+include <tbset.h>
+
+#* HISTORY *
+#* B.Simon 24-Aug-1994 original
+# Phil Hodge 8-Apr-1999 Call tbfpri.
+
+# TMATCH -- Find closest matching rows between two tables
+
+procedure tmatch ()
+
+#--
+pointer input1 # First input table
+pointer input2 # Second input table
+pointer output # Output table
+pointer match1 # Columns from first table used to match
+pointer match2 # Columns from second table used to match
+double maxnorm # Maximum value of norm for allowed match
+pointer incol1 # Columns from first table copied to output
+pointer incol2 # Columns from second table copied to output
+pointer factor # Multiplicative factors used in computing norm
+pointer diagfile # Diagnostic output file
+pointer nmcol1 # Columns from first table in diagnostic output
+pointer nmcol2 # Columns from second table in diagnostic output
+bool sphere # Apply spherical correction to first column?
+
+bool fold
+int mxcol1, mxcol2, ncol1, ncol2, nrow1, nrow2
+int phu_copied # set by tbfpri and ignored
+pointer sp, in1, in2, col1, col2, index1, index2, weight, dist, closest
+
+data fold / false /
+
+string mismatch "Both lists of match columns must have same length"
+string nomatch "Match columns not found in table"
+
+bool clgetb()
+double clgetd()
+int tbpsta()
+pointer tbtopn()
+
+begin
+ # Allocate memory for strings
+
+ call smark (sp)
+ call salloc (input1, SZ_FNAME, TY_CHAR)
+ call salloc (input2, SZ_FNAME, TY_CHAR)
+ call salloc (output, SZ_FNAME, TY_CHAR)
+ call salloc (match1, SZ_FNAME, TY_CHAR)
+ call salloc (match2, SZ_FNAME, TY_CHAR)
+ call salloc (incol1, SZ_FNAME, TY_CHAR)
+ call salloc (incol2, SZ_FNAME, TY_CHAR)
+ call salloc (factor, SZ_FNAME, TY_CHAR)
+ call salloc (diagfile, SZ_FNAME, TY_CHAR)
+ call salloc (nmcol1, SZ_FNAME, TY_CHAR)
+ call salloc (nmcol2, SZ_FNAME, TY_CHAR)
+
+ # Read task parameters
+
+ call clgstr ("input1", Memc[input1], SZ_FNAME)
+ call clgstr ("input2", Memc[input2], SZ_FNAME)
+ call clgstr ("output", Memc[output], SZ_FNAME)
+ call clgstr ("match1", Memc[match1], SZ_FNAME)
+ call clgstr ("match2", Memc[match2], SZ_FNAME)
+ maxnorm = clgetd ("maxnorm")
+
+ call clgstr ("incol1", Memc[incol1], SZ_FNAME)
+ call clgstr ("incol2", Memc[incol2], SZ_FNAME)
+ call clgstr ("factor", Memc[factor], SZ_FNAME)
+ call clgstr ("diagfile", Memc[diagfile], SZ_FNAME)
+ call clgstr ("nmcol1", Memc[nmcol1], SZ_FNAME)
+ call clgstr ("nmcol2", Memc[nmcol2], SZ_FNAME)
+ sphere = clgetb ("sphere")
+
+ # Open input tables and get list of match colums
+
+ in1 = tbtopn (Memc[input1], READ_ONLY, NULL)
+ in2 = tbtopn (Memc[input2], READ_ONLY, NULL)
+
+ mxcol1 = tbpsta (in1, TBL_NCOLS)
+ mxcol2 = tbpsta (in2, TBL_NCOLS)
+
+ call salloc (col1, mxcol1, TY_INT)
+ call salloc (col2, mxcol2, TY_INT)
+
+ call tctexp (in1, Memc[match1], mxcol1, ncol1, Memi[col1])
+ call tctexp (in2, Memc[match2], mxcol2, ncol2, Memi[col2])
+
+ if (ncol1 != ncol2)
+ call error (1, mismatch)
+
+ if (ncol1 == 0)
+ call error (1, nomatch)
+
+ if (ncol1 < 2)
+ sphere = false
+
+ # Sort input tables
+
+ call allrows (in1, nrow1, index1)
+ call allrows (in2, nrow2, index2)
+
+ call tbtsrt (in1, ncol1, Memi[col1], fold, nrow1, Memi[index1])
+ call tbtsrt (in2, ncol2, Memi[col2], fold, nrow2, Memi[index2])
+
+ call salloc (weight, ncol1, TY_DOUBLE)
+ call salloc (dist, nrow1, TY_DOUBLE)
+ call salloc (closest, nrow1, TY_INT)
+
+ # Compute weights from list of factors or table column units
+
+ call getweight (ncol1, Memi[col1], Memi[col2],
+ Memc[factor], Memd[weight])
+
+ # Compute closest match between the two tables
+
+ call getmatch (in1, in2, ncol1, Memi[col1], Memi[col2], Memd[weight],
+ nrow1, Memi[index1], nrow2, Memi[index2], maxnorm,
+ sphere, Memi[closest], Memd[dist])
+
+ # Write output table
+
+ call tbfpri (Memc[input1], Memc[output], phu_copied)
+ call putmatch (Memc[output], Memc[incol1], Memc[incol2], in1, in2,
+ nrow1, Memi[closest])
+
+ # Write diagnostic info
+
+ call infomatch (Memc[diagfile], in1, in2, Memc[nmcol1], Memc[nmcol2],
+ maxnorm, nrow1, Memi[closest], Memd[dist])
+
+ # Clean up
+
+ call mfree (index1, TY_INT)
+ call mfree (index2, TY_INT)
+ call tbtclo (in1)
+ call tbtclo (in2)
+ call sfree (sp)
+end
generated by cgit (git 2.34.1) at 2025-09-20 08:44:11 -0400