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
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
|
include <imhdr.h>
include <error.h>
define XYZ 1 # xyz -> xyz (identity)
define XZY 2 # xyz -> xzy
define YXZ 3 # xyz -> yxz
define YZX 4 # xyz -> yzx
define ZXY 5 # xyz -> zxy
define ZYX 6 # xyz -> zyx
# T_IM3DTRAN -- Transpose 3d images.
#
# The input and output images are given by image template lists. The
# number of output images must match the number of input images. Image
# sections are allowed in the input images and are ignored in the output
# images. If the input and output image names are the same then the transpose
# is performed to a temporary file which then replaces the input image.
procedure t_im3dtran ()
char imtlist1[SZ_LINE] # Input image list
char imtlist2[SZ_LINE] # Output image list
int len_blk # 1D length of transpose block
char image1[SZ_FNAME] # Input image name
char image2[SZ_FNAME] # Output image name
char imtemp[SZ_FNAME] # Temporary file
int list1, list2, new_ax[3], which3d
pointer im1, im2
int clgeti(), imtopen(), imtgetim(), imtlen(), whichtran()
pointer immap()
begin
# Get input and output image template lists, the size of the transpose
# block, and the transpose mapping.
call clgstr ("input", imtlist1, SZ_LINE)
call clgstr ("output", imtlist2, SZ_LINE)
len_blk = clgeti ("len_blk")
new_ax[1] = clgeti ("new_x")
new_ax[2] = clgeti ("new_y")
new_ax[3] = clgeti ("new_z")
# Determine the type of 3d transpose.
which3d = whichtran (new_ax)
if (which3d <= 0)
call error (0, "Invalid mapping of new_x, new_y, new_z")
# Expand the input and output image lists.
list1 = imtopen (imtlist1)
list2 = imtopen (imtlist2)
if (imtlen (list1) != imtlen (list2)) {
call imtclose (list1)
call imtclose (list2)
call error (1, "Number of input and output images not the same")
}
# Do each set of input/output images.
while ((imtgetim (list1, image1, SZ_FNAME) != EOF) &&
(imtgetim (list2, image2, SZ_FNAME) != EOF)) {
call xt_mkimtemp (image1, image2, imtemp, SZ_FNAME)
im1 = immap (image1, READ_ONLY, 0)
im2 = immap (image2, NEW_COPY, im1)
# Do the transpose.
call im3dtranspose (im1, im2, len_blk, which3d, new_ax)
# Unmap the input and output images.
call imunmap (im1)
call imunmap (im2)
call xt_delimtemp (image2, imtemp)
}
call imtclose (list1)
call imtclose (list2)
end
# IM3DTRANSPOSE -- Transpose an image.
#
# Divide the image into square blocks of size len_blk by len_blk.
# Transpose each block with a generic array transpose operator.
procedure im3dtranspose (im_in, im_out, len_blk, which3d, new_ax)
pointer im_in # Input image descriptor
pointer im_out # Output image descriptor
int len_blk # 1D length of transpose block
int which3d # Parameterized transpose order
int new_ax[3] # Map old axis[index] to new value
int x1, x2, nx
int y1, y2, ny
int z1, z2, nz
pointer buf_in, buf_out
pointer imgs3s(), imps3s(), imgs3i(), imps3i(), imgs3l(), imps3l()
pointer imgs3r(), imps3r(), imgs3d(), imps3d(), imgs3x(), imps3x()
begin
# Output image is a copy of input image with dims transposed.
IM_LEN (im_out, 1) = IM_LEN (im_in, new_ax[1])
IM_LEN (im_out, 2) = IM_LEN (im_in, new_ax[2])
IM_LEN (im_out, 3) = IM_LEN (im_in, new_ax[3])
# Break the input image into blocks of at most (len_blk)**3 .
do x1 = 1, IM_LEN (im_in, 1), len_blk {
x2 = x1 + len_blk - 1
if (x2 > IM_LEN(im_in, 1))
x2 = IM_LEN(im_in, 1)
nx = x2 - x1 + 1
do y1 = 1, IM_LEN (im_in, 2), len_blk {
y2 = y1 + len_blk - 1
if (y2 > IM_LEN(im_in, 2))
y2 = IM_LEN(im_in, 2)
ny = y2 - y1 + 1
do z1 = 1, IM_LEN (im_in, 3), len_blk {
z2 = z1 + len_blk - 1
if (z2 > IM_LEN(im_in, 3))
z2 = IM_LEN(im_in, 3)
nz = z2 - z1 + 1
# Switch on the pixel type to optimize IMIO.
switch (IM_PIXTYPE (im_in)) {
case TY_SHORT:
buf_in = imgs3s (im_in, x1, x2, y1, y2, z1, z2)
switch (which3d) {
case XYZ:
buf_out = imps3s (im_out, x1, x2, y1, y2, z1, z2)
call txyz3s (Mems[buf_in], Mems[buf_out], nx,ny,nz)
case XZY:
buf_out = imps3s (im_out, x1, x2, z1, z2, y1, y2)
call txzy3s (Mems[buf_in], Mems[buf_out], nx,ny,nz)
case YXZ:
buf_out = imps3s (im_out, y1, y2, x1, x2, z1, z2)
call tyxz3s (Mems[buf_in], Mems[buf_out], nx,ny,nz)
case YZX:
buf_out = imps3s (im_out, y1, y2, z1, z2, x1, x2)
call tyzx3s (Mems[buf_in], Mems[buf_out], nx,ny,nz)
case ZXY:
buf_out = imps3s (im_out, z1, z2, x1, x2, y1, y2)
call tzxy3s (Mems[buf_in], Mems[buf_out], nx,ny,nz)
case ZYX:
buf_out = imps3s (im_out, z1, z2, y1, y2, x1, x2)
call tzyx3s (Mems[buf_in], Mems[buf_out], nx,ny,nz)
}
case TY_INT:
buf_in = imgs3i (im_in, x1, x2, y1, y2, z1, z2)
switch (which3d) {
case XYZ:
buf_out = imps3i (im_out, x1, x2, y1, y2, z1, z2)
call txyz3i (Memi[buf_in], Memi[buf_out], nx,ny,nz)
case XZY:
buf_out = imps3i (im_out, x1, x2, z1, z2, y1, y2)
call txzy3i (Memi[buf_in], Memi[buf_out], nx,ny,nz)
case YXZ:
buf_out = imps3i (im_out, y1, y2, x1, x2, z1, z2)
call tyxz3i (Memi[buf_in], Memi[buf_out], nx,ny,nz)
case YZX:
buf_out = imps3i (im_out, y1, y2, z1, z2, x1, x2)
call tyzx3i (Memi[buf_in], Memi[buf_out], nx,ny,nz)
case ZXY:
buf_out = imps3i (im_out, z1, z2, x1, x2, y1, y2)
call tzxy3i (Memi[buf_in], Memi[buf_out], nx,ny,nz)
case ZYX:
buf_out = imps3i (im_out, z1, z2, y1, y2, x1, x2)
call tzyx3i (Memi[buf_in], Memi[buf_out], nx,ny,nz)
}
case TY_LONG:
buf_in = imgs3l (im_in, x1, x2, y1, y2, z1, z2)
switch (which3d) {
case XYZ:
buf_out = imps3l (im_out, x1, x2, y1, y2, z1, z2)
call txyz3l (Meml[buf_in], Meml[buf_out], nx,ny,nz)
case XZY:
buf_out = imps3l (im_out, x1, x2, z1, z2, y1, y2)
call txzy3l (Meml[buf_in], Meml[buf_out], nx,ny,nz)
case YXZ:
buf_out = imps3l (im_out, y1, y2, x1, x2, z1, z2)
call tyxz3l (Meml[buf_in], Meml[buf_out], nx,ny,nz)
case YZX:
buf_out = imps3l (im_out, y1, y2, z1, z2, x1, x2)
call tyzx3l (Meml[buf_in], Meml[buf_out], nx,ny,nz)
case ZXY:
buf_out = imps3l (im_out, z1, z2, x1, x2, y1, y2)
call tzxy3l (Meml[buf_in], Meml[buf_out], nx,ny,nz)
case ZYX:
buf_out = imps3l (im_out, z1, z2, y1, y2, x1, x2)
call tzyx3l (Meml[buf_in], Meml[buf_out], nx,ny,nz)
}
case TY_REAL:
buf_in = imgs3r (im_in, x1, x2, y1, y2, z1, z2)
switch (which3d) {
case XYZ:
buf_out = imps3r (im_out, x1, x2, y1, y2, z1, z2)
call txyz3r (Memr[buf_in], Memr[buf_out], nx,ny,nz)
case XZY:
buf_out = imps3r (im_out, x1, x2, z1, z2, y1, y2)
call txzy3r (Memr[buf_in], Memr[buf_out], nx,ny,nz)
case YXZ:
buf_out = imps3r (im_out, y1, y2, x1, x2, z1, z2)
call tyxz3r (Memr[buf_in], Memr[buf_out], nx,ny,nz)
case YZX:
buf_out = imps3r (im_out, y1, y2, z1, z2, x1, x2)
call tyzx3r (Memr[buf_in], Memr[buf_out], nx,ny,nz)
case ZXY:
buf_out = imps3r (im_out, z1, z2, x1, x2, y1, y2)
call tzxy3r (Memr[buf_in], Memr[buf_out], nx,ny,nz)
case ZYX:
buf_out = imps3r (im_out, z1, z2, y1, y2, x1, x2)
call tzyx3r (Memr[buf_in], Memr[buf_out], nx,ny,nz)
}
case TY_DOUBLE:
buf_in = imgs3d (im_in, x1, x2, y1, y2, z1, z2)
switch (which3d) {
case XYZ:
buf_out = imps3d (im_out, x1, x2, y1, y2, z1, z2)
call txyz3d (Memd[buf_in], Memd[buf_out], nx,ny,nz)
case XZY:
buf_out = imps3d (im_out, x1, x2, z1, z2, y1, y2)
call txzy3d (Memd[buf_in], Memd[buf_out], nx,ny,nz)
case YXZ:
buf_out = imps3d (im_out, y1, y2, x1, x2, z1, z2)
call tyxz3d (Memd[buf_in], Memd[buf_out], nx,ny,nz)
case YZX:
buf_out = imps3d (im_out, y1, y2, z1, z2, x1, x2)
call tyzx3d (Memd[buf_in], Memd[buf_out], nx,ny,nz)
case ZXY:
buf_out = imps3d (im_out, z1, z2, x1, x2, y1, y2)
call tzxy3d (Memd[buf_in], Memd[buf_out], nx,ny,nz)
case ZYX:
buf_out = imps3d (im_out, z1, z2, y1, y2, x1, x2)
call tzyx3d (Memd[buf_in], Memd[buf_out], nx,ny,nz)
}
case TY_COMPLEX:
buf_in = imgs3x (im_in, x1, x2, y1, y2, z1, z2)
switch (which3d) {
case XYZ:
buf_out = imps3x (im_out, x1, x2, y1, y2, z1, z2)
call txyz3x (Memx[buf_in], Memx[buf_out], nx,ny,nz)
case XZY:
buf_out = imps3x (im_out, x1, x2, z1, z2, y1, y2)
call txzy3x (Memx[buf_in], Memx[buf_out], nx,ny,nz)
case YXZ:
buf_out = imps3x (im_out, y1, y2, x1, x2, z1, z2)
call tyxz3x (Memx[buf_in], Memx[buf_out], nx,ny,nz)
case YZX:
buf_out = imps3x (im_out, y1, y2, z1, z2, x1, x2)
call tyzx3x (Memx[buf_in], Memx[buf_out], nx,ny,nz)
case ZXY:
buf_out = imps3x (im_out, z1, z2, x1, x2, y1, y2)
call tzxy3x (Memx[buf_in], Memx[buf_out], nx,ny,nz)
case ZYX:
buf_out = imps3x (im_out, z1, z2, y1, y2, x1, x2)
call tzyx3x (Memx[buf_in], Memx[buf_out], nx,ny,nz)
}
default:
call error (3, "unknown pixel type")
}
}
}
}
end
# WHICHTRAN -- Return transpose type.
int procedure whichtran (new_ax)
int new_ax[3]
int which
begin
which = 0
if (new_ax[1] == 1) {
if (new_ax[2] == 2)
which = XYZ
else if (new_ax[2] == 3)
which = XZY
} else if (new_ax[1] == 2) {
if (new_ax[2] == 1)
which = YXZ
else if (new_ax[2] == 3)
which = YZX
} else if (new_ax[1] == 3) {
if (new_ax[2] == 1)
which = ZXY
else if (new_ax[2] == 2)
which = ZYX
}
return (which)
end
|