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
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
|
include <imhdr.h>
include "quadgeom.h"
define OPPERATIONS "|multiply|divide|"
define OPMULTIPLY 1
define OPDIVIDE 2
define TOL1 1.0e-7
define TOL2 -1.0e-7
procedure t_quadscale ()
char input[SZ_FNAME] #TI Input image name.
char output[SZ_FNAME] #TI Output image name.
char instrument[SZ_FNAME] #TI Instrument translation file
real gain[QG_MAXAMPS] #TI Gain factor for each quadrant
int op #TI Multiply or divide by gain factors
int i
pointer in, out, qg
char buffer[SZ_LINE]
real clgetr()
int clgwrd(), hdmaccf()
pointer immap()
begin
# Open instrument file
call clgstr ("instrument", instrument, SZ_FNAME)
call hdmopen (instrument)
# Open input image
call clgstr ("input", input, SZ_FNAME)
in = immap (input, READ_ONLY, 0)
# Allocate quadgeom structure and initialise it from image header
call quadalloc (qg)
# Fill out quadgeom structure from header depending on revision level
if (hdmaccf (in, "HDR_REV") == NO) {
# Check to see if the image has been processed or not
if (hdmaccf (in, "ccdproc") == YES) {
call quadgeomred (in, qg)
} else {
call quadgeom (in, qg, "", "")
}
} else {
call qghdr2 (in, qg)
}
# Open output image
call clgstr ("output", output, SZ_FNAME)
out = immap (output, NEW_COPY, in)
IM_PIXTYPE(out) = TY_REAL
# Get gain factors
gain[1] = clgetr ("gain11")
gain[2] = clgetr ("gain12")
gain[3] = clgetr ("gain21")
gain[4] = clgetr ("gain22")
# Get direction of opperation
op = clgwrd ("opperation", buffer, SZ_LINE, OPPERATIONS)
switch (op) {
case OPMULTIPLY:
call quadmult (in, out, gain, qg)
case OPDIVIDE:
# Check for zero gain --> divide by zero
do i = 1, 4 {
if ((gain[i] < TOL1) && (gain[i] > TOL2)) {
call error (0, "Attempt to divide by gain value of zero")
}
}
call quaddiv (in, out, gain, qg)
}
call imunmap (in)
call imunmap (out)
call quadfree (qg)
call hdmclose ()
end
procedure quadmult (in, out, gain, qg)
pointer in #I imio pointer for input image.
pointer out #I imio pointer for output image.
real gain[ARB] #I Array of gain factors.
pointer qg #I Pointer to quadgeom structure.
pointer inbuf, obuf
int junk, nx, x, y, line, amp, amp2, off
long invec[IM_MAXDIM], ovec[IM_MAXDIM]
int imgnlr(), impnlr()
begin
# Setup start vector for sequential reads
call amovkl (long(1), invec, IM_MAXDIM)
call amovkl (long(1), ovec, IM_MAXDIM)
do y = 1, QG_NAMPS(qg) {
amp2 = QG_AMP(qg, 1, y)
do line = 1, QG_NY(qg, amp2) {
junk = imgnlr (in, inbuf, invec)
junk = impnlr (out, obuf, ovec)
off = 0
do x = 1, QG_NAMPSX(qg) {
amp = QG_AMP(qg, x, y)
nx = QG_NX(qg, amp)
call amulkr (Memr[inbuf+off], gain[amp], Memr[obuf+off], nx)
off = off + nx
}
}
}
end
procedure quaddiv (in, out, gain, qg)
pointer in #I imio pointer for input image.
pointer out #I imio pointer for output image.
real gain[ARB] #I Array of gain factors.
pointer qg #I Pointer to quadgeom structure.
pointer inbuf, obuf
int junk, nx, x, y, line, amp, amp2, off
long invec[IM_MAXDIM], ovec[IM_MAXDIM]
int imgnlr(), impnlr()
begin
# Setup start vector for sequential reads
call amovkl (long(1), invec, IM_MAXDIM)
call amovkl (long(1), ovec, IM_MAXDIM)
do y = 1, QG_NAMPS(qg) {
amp2 = QG_AMP(qg, 1, y)
do line = 1, QG_NY(qg, amp2) {
junk = imgnlr (in, inbuf, invec)
junk = impnlr (out, obuf, ovec)
off = 0
do x = 1, QG_NAMPSX(qg) {
amp = QG_AMP(qg, x, y)
nx = QG_NX(qg, amp)
call adivkr (Memr[inbuf+off], gain[amp], Memr[obuf+off], nx)
off = off + nx
}
}
}
end
|
