1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
|
include <error.h>
include <ctype.h>
include "i2sun.h"
# TR_ULUT -- Generates a look up table from data supplied by user. The
# data is read from a two column text file of intensity, greyscale values.
# The input data are sorted, then mapped to the x range [0-4095]. A
# piecewise linear look up table of 4096 values is then constructed from
# the (x,y) pairs given. A pointer to the look up table, as well as the z1
# and z2 intensity endpoints, is returned.
procedure tr_ulut (fname, z1, z2, lut)
char fname[SZ_FNAME] # Name of file with intensity, greyscale values
real z1 # Intensity mapped to minimum gs value
real z2 # Intensity mapped to maximum gs value
pointer lut # Look up table - pointer is returned
pointer sp, x, y
int nvalues, i, j, x1, x2, y1
real delta_gs, delta_xv, slope
errchk ds_rlut, ds_sort, malloc
begin
call smark (sp)
call salloc (x, U_MAXPTS, TY_REAL)
call salloc (y, U_MAXPTS, TY_REAL)
# Read intensities and greyscales from the user's input file. The
# intensity range is then mapped into a standard range and the
# values sorted.
call ds_rlut (fname, Memr[x], Memr[y], nvalues)
call alimr (Memr[x], nvalues, z1, z2)
call amapr (Memr[x], Memr[x], nvalues, z1, z2, real(U_Z1), real(U_Z2))
call ds_sort (Memr[x], Memr[y], nvalues)
# Fill lut in straight line segments - piecewise linear
call malloc (lut, U_MAXPTS, TY_SHORT)
do i = 1, nvalues-1 {
delta_gs = Memr[y+i] - Memr[y+i-1]
delta_xv = Memr[x+i] - Memr[x+i-1]
slope = delta_gs / delta_xv
x1 = int (Memr[x+i-1])
x2 = int (Memr[x+i])
y1 = int (Memr[y+i-1])
do j = x1, x2
Mems[lut+j] = y1 + slope * (j-x1)
}
Mems[lut+U_MAXPTS-1] = y1 + (slope * U_Z2)
call sfree (sp)
end
# DS_RLUT -- Read text file of x, y, values.
procedure ds_rlut (utab, x, y, nvalues)
char utab[SZ_FNAME] # Name of list file
real x[U_MAXPTS] # Array of x values, filled on return
real y[U_MAXPTS] # Array of y values, filled on return
int nvalues # Number of values in x, y vectors - returned
int n, fd
pointer sp, lbuf, ip
real xval, yval
int getline(), open()
errchk open, sscan, getline, salloc
begin
call smark (sp)
call salloc (lbuf, SZ_LINE, TY_CHAR)
iferr (fd = open (utab, READ_ONLY, TEXT_FILE))
call error (1, "Error opening user lookup table")
n = 0
while (getline (fd, Memc[lbuf]) != EOF) {
# Skip comment lines and blank lines.
if (Memc[lbuf] == '#')
next
for (ip=lbuf; IS_WHITE(Memc[ip]); ip=ip+1)
;
if (Memc[ip] == '\n' || Memc[ip] == EOS)
next
# Decode the points to be plotted.
call sscan (Memc[ip])
call gargr (xval)
call gargr (yval)
n = n + 1
if (n > U_MAXPTS)
call error (2,
"Intensity transformation table cannot exceed 4096 values")
x[n] = xval
y[n] = yval
}
nvalues = n
call close (fd)
call sfree (sp)
end
# DS_SORT -- Bubble sort of paired arrays.
procedure ds_sort (xvals, yvals, nvals)
real xvals[nvals] # Array of x values
real yvals[nvals] # Array of y values
int nvals # Number of values in each array
int i, j
real temp
define swap {temp=$1;$1=$2;$2=temp}
begin
for (i=nvals; i > 1; i=i-1)
for (j=1; j < i; j=j+1)
if (xvals[j] > xvals[j+1]) {
# Out of order; exchange y values
swap (xvals[j], xvals[j+1])
swap (yvals[j], yvals[j+1])
}
end
|