1 /* xpsgip.f -- translated by f2c (version 19980913).
2 You must link the resulting object file with the libraries:
3 -lf2c -lm (in that order)
4 */
5
6 #include "f2c.h"
7
8 /* $Procedure XPSGIP ( Transpose a matrix, general dimension, in place ) */
xpsgip_(integer * nrow,integer * ncol,doublereal * matrix)9 /* Subroutine */ int xpsgip_(integer *nrow, integer *ncol, doublereal *matrix)
10 {
11 integer dest;
12 doublereal temp;
13 integer k, m, n, r__, moved, start;
14 doublereal source;
15 integer nmoves;
16
17 /* $ Abstract */
18
19 /* Transpose a matrix of arbitrary size and shape in place. */
20
21 /* $ Disclaimer */
22
23 /* THIS SOFTWARE AND ANY RELATED MATERIALS WERE CREATED BY THE */
24 /* CALIFORNIA INSTITUTE OF TECHNOLOGY (CALTECH) UNDER A U.S. */
25 /* GOVERNMENT CONTRACT WITH THE NATIONAL AERONAUTICS AND SPACE */
26 /* ADMINISTRATION (NASA). THE SOFTWARE IS TECHNOLOGY AND SOFTWARE */
27 /* PUBLICLY AVAILABLE UNDER U.S. EXPORT LAWS AND IS PROVIDED "AS-IS" */
28 /* TO THE RECIPIENT WITHOUT WARRANTY OF ANY KIND, INCLUDING ANY */
29 /* WARRANTIES OF PERFORMANCE OR MERCHANTABILITY OR FITNESS FOR A */
30 /* PARTICULAR USE OR PURPOSE (AS SET FORTH IN UNITED STATES UCC */
31 /* SECTIONS 2312-2313) OR FOR ANY PURPOSE WHATSOEVER, FOR THE */
32 /* SOFTWARE AND RELATED MATERIALS, HOWEVER USED. */
33
34 /* IN NO EVENT SHALL CALTECH, ITS JET PROPULSION LABORATORY, OR NASA */
35 /* BE LIABLE FOR ANY DAMAGES AND/OR COSTS, INCLUDING, BUT NOT */
36 /* LIMITED TO, INCIDENTAL OR CONSEQUENTIAL DAMAGES OF ANY KIND, */
37 /* INCLUDING ECONOMIC DAMAGE OR INJURY TO PROPERTY AND LOST PROFITS, */
38 /* REGARDLESS OF WHETHER CALTECH, JPL, OR NASA BE ADVISED, HAVE */
39 /* REASON TO KNOW, OR, IN FACT, SHALL KNOW OF THE POSSIBILITY. */
40
41 /* RECIPIENT BEARS ALL RISK RELATING TO QUALITY AND PERFORMANCE OF */
42 /* THE SOFTWARE AND ANY RELATED MATERIALS, AND AGREES TO INDEMNIFY */
43 /* CALTECH AND NASA FOR ALL THIRD-PARTY CLAIMS RESULTING FROM THE */
44 /* ACTIONS OF RECIPIENT IN THE USE OF THE SOFTWARE. */
45
46 /* $ Required_Reading */
47
48 /* None. */
49
50 /* $ Keywords */
51
52 /* MATRIX */
53
54 /* $ Declarations */
55 /* $ Brief_I/O */
56
57 /* VARIABLE I/O DESCRIPTION */
58 /* -------- --- -------------------------------------------------- */
59 /* NROW I Number of rows of input matrix. */
60 /* NCOL I Number of columns of input matrix. */
61 /* MATRIX I-O Matrix to be transposed/transposed matrix. */
62
63 /* $ Detailed_Input */
64
65 /* MATRIX Matrix to be transposed. */
66
67 /* NROW Number of rows of input matrix MATRIX. */
68
69 /* NCOL Number of columns of input matrix MATRIX. */
70
71 /* $ Detailed_Output */
72
73 /* MATRIX Transposed matrix: element (i,j) of the input */
74 /* matrix is element (j,i) of the output matrix. */
75
76 /* $ Parameters */
77
78 /* None. */
79
80 /* $ Exceptions */
81
82 /* Error Free. */
83
84 /* 1) If either NROW or NCOL is less than or equal to zero, no action */
85 /* is taken. The routine simply returns. */
86
87 /* $ Files */
88
89 /* None. */
90
91 /* $ Particulars */
92
93 /* This routine replaces the input matrix MATRIX with its transpose. */
94
95 /* NOTE: The matrix MATRIX is declared one-dimensional for */
96 /* computational purposes only. The calling program may */
97 /* declare it as MATRIX(NROW,NCOL) or MATRIX(NCOL,NROW). */
98
99 /* This routine assumes that the elements of the matrix to be */
100 /* transformed are stored in contiguous memory locations as */
101 /* shown here. On output these elements will be rearranged */
102 /* in consecutive memory locations as shown. */
103
104 /* MATRIX on input MATRIX on output */
105
106 /* m_11 m_11 */
107 /* m_21 m_12 */
108 /* m_31 m_13 */
109 /* . . */
110 /* . . */
111 /* . m_1ncol */
112 /* m_nrow1 m_21 */
113 /* m_12 m_22 */
114 /* m_22 m_23 */
115 /* m_32 . */
116 /* . . */
117 /* . m_2ncol */
118 /* . . */
119 /* m_nrow2 */
120 /* . . */
121
122 /* . . */
123
124 /* . . */
125 /* m_1ncol */
126 /* m_2ncol m_nrow1 */
127 /* m_3ncol m_nrow2 */
128 /* . m_nrow3 */
129 /* . . */
130 /* . . */
131 /* m_nrowncol m_nrowncol */
132
133
134 /* For those familiar with permutations, this algorithm relies */
135 /* upon the fact that the transposition of a matrix, which has */
136 /* been stored as a 1-dimensional array, is simply the action of a */
137 /* permutation applied to that array. Since any permutation */
138 /* can be decomposed as a product of disjoint cycles, it is */
139 /* possible to transpose the matrix with only one additional */
140 /* storage register. However, once a cycle has been computed */
141 /* it is necessary to find the next entry in the array that */
142 /* has not been moved by the permutation. For this reason the */
143 /* algorithm is slower than would be necessary if the numbers */
144 /* of rows and columns were known in advance. */
145
146 /* $ Examples */
147
148 /* This routine is provided for situation where it is convenient to */
149 /* transpose a general two-dimensional matrix */
150 /* in place rather than store the result in a */
151 /* separate array. Note that the call */
152
153 /* CALL XPOSEG ( MATRIX, NROW, NCOL, MATRIX ) */
154
155 /* is not permitted by the ANSI Fortran 77 standard; this routine */
156 /* can be called instead to achieve the same result. */
157
158 /* $ Restrictions */
159
160 /* None. */
161
162 /* $ Author_and_Institution */
163
164 /* N.J. Bachman (JPL) */
165 /* W.L. Taber (JPL) */
166
167 /* $ Literature_References */
168
169 /* None. */
170
171 /* $ Version */
172
173 /* - SPICELIB Version 1.0.1, 19-SEP-2006 (EDW) */
174
175 /* Initial version date unknown. Version data entry */
176 /* added this date. */
177
178
179 /* -& */
180 /* $ Index_Entries */
181
182 /* transpose a matrix general */
183
184 /* -& */
185
186 /* Local Variables */
187
188
189 /* Take care of dumb cases first. */
190
191 if (*nrow <= 0 || *ncol <= 0) {
192 return 0;
193 }
194
195 /* Use the abbreviations N and M for NROW and NCOL. */
196
197 n = *nrow;
198 m = *ncol;
199
200 /* Set up the upper bound for the number of objects to be moved and */
201 /* initialize the counters. */
202
203 nmoves = n * m - 2;
204 moved = 0;
205 start = 1;
206
207 /* Until MOVED is equal to NMOVES, there is some matrix element that */
208 /* has not been moved to its proper location in the transpose matrix. */
209
210 while(moved < nmoves) {
211 source = matrix[start];
212 k = start / n;
213 r__ = start - n * k;
214 dest = r__ * m + k;
215
216 /* Perform this cycle of the permutation. We will be done when */
217 /* the destination of the next element is equal to the starting */
218 /* position of the first element to be moved in this cycle. */
219
220 while(dest != start) {
221 temp = matrix[dest];
222 matrix[dest] = source;
223 source = temp;
224 ++moved;
225 k = dest / n;
226 r__ = dest - k * n;
227 dest = m * r__ + k;
228 }
229 matrix[dest] = source;
230 dest = 0;
231 ++moved;
232
233 /* Find the next element of the matrix that has not already been */
234 /* moved by the transposition operation. */
235
236 if (moved < nmoves) {
237 while(dest != start) {
238 ++start;
239 k = start / n;
240 r__ = start - k * n;
241 dest = r__ * m + k;
242 while(dest > start) {
243 k = dest / n;
244 r__ = dest - k * n;
245 dest = m * r__ + k;
246 }
247 }
248 }
249 }
250 return 0;
251 } /* xpsgip_ */
252
253