1 #include "mrilib.h"
2
3 /** NOT 7D SAFE **/
4
5 /*----------------------------------------------------------------------*/
6 /*! Flip a 3D image. The (outx,outy,outz) parameters specify the
7 direction of the output axes relative to the input axes:
8 - +1 => +x
9 - -1 => -x
10 - +2 => +y
11 - -2 => -y
12 - +3 => +z
13 - -3 => -z
14
15 Among (outx,outy,outz), exactly one must be either +1 or -1, one
16 must be +2 or -2, and one must be +3 or -3. Bad inputs result in
17 a NULL return value.
18
19 If the input image doesn't have a data array, then the output
20 image won't either.
21 ------------------------------------------------------------------------*/
22
mri_flip3D(int outx,int outy,int outz,MRI_IMAGE * inim)23 MRI_IMAGE * mri_flip3D( int outx, int outy, int outz, MRI_IMAGE *inim )
24 {
25 MRI_IMAGE *outim ;
26 int ii,jj,kk , dsiz , nxin,nyin,nzin , nxout,nyout,nzout ;
27 int nxyin , nxyout , ijk_out , ijk_in ;
28 int ax,bx,cx,dx , ay,by,cy,dy , az,bz,cz,dz , aa,bb,cc,dd ;
29 char *inar , *outar ;
30 float delx,dely,delz ;
31
32 ENTRY("mri_flip3D") ;
33
34 /* check inputs for correctness */
35
36 if( inim == NULL || outx == 0 || outy == 0 || outz == 0 ) RETURN( NULL );
37 ii = abs(outx) ; jj = abs(outy) ; kk = abs(outz) ;
38 if( ii > 3 || jj > 3 || kk > 3 ) RETURN( NULL );
39 if( ii == jj || ii == kk || jj == kk ) RETURN( NULL );
40 if( ii+jj+kk != 6 ) RETURN( NULL );
41
42 if( outx==1 && outy==2 && outz==3 ) RETURN( mri_copy(inim) ); /* easy case */
43
44 nxin = inim->nx ;
45 nyin = inim->ny ; nxyin = nxin*nyin ;
46 nzin = inim->nz ;
47
48 /* setup so that i_out = ax + bx*i_in + cx*j_in + dx*k_in,
49 for i_in=0..nxin-1, j_in=0..nyin-1, k_in=0..nzin-1,
50 and then similarly for y and z axes */
51
52 switch( outx ){
53 case 1: ax=0 ; bx= 1; cx= 0; dx= 0; nxout=nxin; delx=inim->dx; break;
54 case -1: ax=nxin-1; bx=-1; cx= 0; dx= 0; nxout=nxin; delx=inim->dx; break;
55 case 2: ax=0 ; bx= 0; cx= 1; dx= 0; nxout=nyin; delx=inim->dy; break;
56 case -2: ax=nyin-1; bx= 0; cx=-1; dx= 0; nxout=nyin; delx=inim->dy; break;
57 case 3: ax=0 ; bx= 0; cx= 0; dx= 1; nxout=nzin; delx=inim->dz; break;
58 case -3: ax=nzin-1; bx= 0; cx= 0; dx=-1; nxout=nzin; delx=inim->dz; break;
59 default: RETURN( NULL );
60 }
61 switch( outy ){
62 case 1: ay=0 ; by= 1; cy= 0; dy= 0; nyout=nxin; dely=inim->dx; break;
63 case -1: ay=nxin-1; by=-1; cy= 0; dy= 0; nyout=nxin; dely=inim->dx; break;
64 case 2: ay=0 ; by= 0; cy= 1; dy= 0; nyout=nyin; dely=inim->dy; break;
65 case -2: ay=nyin-1; by= 0; cy=-1; dy= 0; nyout=nyin; dely=inim->dy; break;
66 case 3: ay=0 ; by= 0; cy= 0; dy= 1; nyout=nzin; dely=inim->dz; break;
67 case -3: ay=nzin-1; by= 0; cy= 0; dy=-1; nyout=nzin; dely=inim->dz; break;
68 default: RETURN( NULL );
69 }
70 switch( outz ){
71 case 1: az=0 ; bz= 1; cz= 0; dz= 0; nzout=nxin; delz=inim->dx; break;
72 case -1: az=nxin-1; bz=-1; cz= 0; dz= 0; nzout=nxin; delz=inim->dx; break;
73 case 2: az=0 ; bz= 0; cz= 1; dz= 0; nzout=nyin; delz=inim->dy; break;
74 case -2: az=nyin-1; bz= 0; cz=-1; dz= 0; nzout=nyin; delz=inim->dy; break;
75 case 3: az=0 ; bz= 0; cz= 0; dz= 1; nzout=nzin; delz=inim->dz; break;
76 case -3: az=nzin-1; bz= 0; cz= 0; dz=-1; nzout=nzin; delz=inim->dz; break;
77 default: RETURN( NULL );
78 }
79 nxyout = nxout*nyout ;
80
81 /* 3D index ijk_out = i_out + nxout*j_out + nxyout*k_out
82
83 = (ax + bx*i_in + cx*j_in + dx*k_in)
84 +(ay + by*i_in + cy*j_in + dy*k_in)*nxout
85 +(az + bz*i_in + cz*j_in + dz*k_in)*nxyout
86
87 = (ax+ay*nxout+az*nxyout)
88 +(bx+by*nxout+bz*nxyout)*i_in
89 +(cx+cy*nxout+cz*nxyout)*j_in
90 +(dx+dy*nxout+dz*nxyout)*k_in
91
92 = aa + bb*i_in + cc*j_in + dd*k_in */
93
94 inar = mri_data_pointer( inim ) ;
95 outim = mri_new_vol( nxout,nyout,nzout , inim->kind ) ;
96
97 outim->dx = delx ; outim->dy = dely ; outim->dz = delz ;
98
99 outar = mri_data_pointer( outim ) ;
100
101 if( inar == NULL ){ /* empty input ==> empty output */
102 free(outar) ; mri_fix_data_pointer(NULL,outim ); RETURN(outim);
103 }
104
105 dsiz = outim->pixel_size ; /* size of each voxel in bytes */
106
107 aa = (ax+ay*nxout+az*nxyout)*dsiz ; /* same as aa, etc. in above */
108 bb = (bx+by*nxout+bz*nxyout)*dsiz ; /* comment, but scaled to be */
109 cc = (cx+cy*nxout+cz*nxyout)*dsiz ; /* address offset in bytes */
110 dd = (dx+dy*nxout+dz*nxyout)*dsiz ;
111
112 for( ijk_in=kk=0 ; kk < nzin ; kk++ ){
113 for( jj=0 ; jj < nyin ; jj++ ){
114 for( ii=0 ; ii < nxin ; ii++,ijk_in+=dsiz ){
115 ijk_out = aa + bb*ii + cc*jj + dd*kk ;
116 memcpy( outar+ijk_out , inar+ijk_in , dsiz ) ;
117 }
118 }
119 }
120
121 RETURN(outim) ;
122 }
123