~ [ source navigation ] ~ [ diff markup ] ~ [ identifier search ] ~ [ freetext search ] ~ [ file search ] ~

Linux Cross Reference
Tina6/tina-libs/tina/math/mathMatv_mat.c

Version: ~
Architecture: ~ [ i386 ] ~ [ alpha ] ~ [ m68k ] ~ [ mips ] ~ [ ppc ] ~ [ sparc ] ~ [ sparc64 ] ~

  1 /**********
  2  *
  3  * Copyright (c) 2003, Division of Imaging Science and Biomedical Engineering,
  4  * University of Manchester, UK.  All rights reserved.
  5  * 
  6  * Redistribution and use in source and binary forms, with or without modification, 
  7  * are permitted provided that the following conditions are met:
  8  * 
  9  *   . Redistributions of source code must retain the above copyright notice, 
 10  *     this list of conditions and the following disclaimer.
 11  *    
 12  *   . Redistributions in binary form must reproduce the above copyright notice,
 13  *     this list of conditions and the following disclaimer in the documentation 
 14  *     and/or other materials provided with the distribution.
 15  * 
 16  *   . Neither the name of the University of Manchester nor the names of its
 17  *     contributors may be used to endorse or promote products derived from this 
 18  *     software without specific prior written permission.
 19  * 
 20  * 
 21  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" 
 22  * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 
 23  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 
 24  * ARE DISCLAIMED.  IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE 
 25  * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR 
 26  * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF 
 27  * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
 28  * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN 
 29  * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) 
 30  * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE 
 31  * POSSIBILITY OF SUCH DAMAGE.
 32  *
 33  **********
 34  *
 35  * Program :    TINA
 36  * File    :  $Source: /home/tina/cvs/tina-libs/tina/math/mathMatv_mat.c,v $
 37  * Date    :  $Date: 2008/12/04 15:28:02 $
 38  * Version :  $Revision: 1.6 $
 39  * CVS Id  :  $Id: mathMatv_mat.c,v 1.6 2008/12/04 15:28:02 paul Exp $
 40  *
 41  * Author  :  Legacy TINA
 42  *
 43  * Notes : Dynamically allocated matrix of doubles.
 44  *
 45  *********
 46 */
 47 /** 
 48  *  @file
 49  *  @brief  Functions pertaining to Mat's, dynamically allocated matrices of doubles.   
 50  *
 51  *  Functions include;
 52  *
 53  *  - Mat allocation and initialisation.
 54  *  - Calculations using Mat's (eg, Mat addition, multiplication).
 55  *  - I/O functions.
 56  *  - Mat access/set/get functions.
 57  *  - conversion to/from Mat3.
 58 */
 59 
 60 #include "mathMatv_mat.h"
 61 
 62 #if HAVE_CONFIG_H
 63 #include <config.h>
 64 #endif
 65 
 66 #include <tina/sys/sysDef.h>
 67 #include <tina/sys/sysPro.h>
 68 #include <tina/math/math_MatvDef.h>
 69 #include <tina/math/math_GeomDef.h>
 70 #include <tina/math/math_GeomPro.h>
 71 #include <tina/math/math_UtilPro.h>
 72 
 73 Mat            *mat_make(int m, int n)
 74 {
 75         Mat            *a = talloc(Mat);
 76         a->m = m;
 77         a->n = n;
 78         a->el = tarray_alloc(0, 0, m, n, double);
 79         return (a);
 80 }
 81 
 82 void            mat_free(Mat * a)
 83 {
 84         if (a == NULL)
 85                 return;
 86         tarray_free(a->el, 0, 0, a->m, a->n, double);
 87         rfree(a);
 88 }
 89 
 90 void            mat_copy(Mat * a, Mat * b)
 91 {
 92         if (a == NULL)
 93                 return;
 94         tarray_copy_inplace(a->el, b->el, 0, 0, a->m, a->n, double);
 95 }
 96 
 97 void            mat_zero(Mat * a)
 98 {
 99         int             i, j;
100         if (a == NULL)
101                 return;
102         for (i = 0; i < a->m; i++)
103                 for (j = 0; j < a->n; j++)
104                         a->el[i][j] = 0;
105 }
106 
107 void            mat_unit(Mat * a)
108 {
109         int             i, j;
110         if (a == NULL)
111                 return;
112         for (i = 0; i < a->m; i++)
113                 for (j = 0; j < a->n; j++)
114                         a->el[i][j] = (i == j) ? 1 : 0;
115 }
116 
117 void            mat_rand_unif(Mat * a, double p, double q)
118 {
119         int             i, j;
120         if (a == NULL)
121                 return;
122         for (i = 0; i < a->m; i++)
123                 for (j = 0; j < a->n; j++)
124                         a->el[i][j] = rand_unif(p, q);
125 }
126 
127 void            mat_rand_normal(Mat * a, double m, double s)
128 {
129         int             i, j;
130         if (a == NULL)
131                 return;
132         for (i = 0; i < a->m; i++)
133                 for (j = 0; j < a->n; j++)
134                         a->el[i][j] = rand_normal(m, s);
135 }
136 
137 void            mat_transp(Mat * at, Mat * a)
138 {
139         int             i, j;
140         if (at == NULL || a == NULL)
141                 return;
142         for (i = 0; i < a->m; i++)
143                 for (j = 0; j < a->n; j++)
144                         at->el[j][i] = a->el[i][j];
145 }
146 
147 void            mat_stransp(Mat * a)
148 {
149         int             i, j;
150         if (a == NULL)
151                 return;
152         for (i = 0; i < a->m; i++)
153                 for (j = 0; j < i; j++)
154                         SWAP(double, a->el[i][j], a->el[j][i]);
155 }
156 
157 void            mat_prod(Mat * ab, Mat * a, Mat * b)
158 {
159         int             i, j, k;
160         if (ab == NULL || a == NULL || b == NULL)
161                 return;
162         for (i = 0; i < a->m; i++)
163                 for (j = 0; j < b->n; j++)
164                 {
165                         double          sum = 0;
166                         for (k = 0; k < a->n; k++)
167                                 sum += a->el[i][k] * b->el[k][j];
168                         ab->el[i][j] = sum;
169                 }
170 }
171 
172 void            mat_sqr(Mat * aa, Mat * a)
173 {
174         int             i, j, k;
175         if (aa == NULL || a == NULL)
176                 return;
177         for (i = 0; i < a->m; i++)
178                 for (j = 0; j < a->m; j++)
179                 {
180                         double          sum = 0;
181                         for (k = 0; k < a->n; k++)
182                                 sum += a->el[i][k] * a->el[j][k];
183                         aa->el[i][j] = sum;
184                 }
185 }
186 
187 void            mat_vprod(Vec * av, Mat * a, Vec * v)
188 {
189         int             i, j;
190         if (av == NULL || a == NULL || v == NULL)
191                 return;
192         for (i = 0; i < a->m; i++)
193         {
194                 double          sum = 0;
195                 for (j = 0; j < a->n; j++)
196                         sum += a->el[i][j] * v->el[j];
197                 av->el[i] = sum;
198         }
199 }
200 
201 void            mat_dprod(Vec * d, Mat * a)
202 {
203         int             i, j;
204         if (a == NULL || d == NULL)
205                 return;
206         for (i = 0; i < a->m; i++)
207                 for (j = 0; j < a->n; j++)
208                         a->el[i][j] *= d->el[i];
209 }
210 
211 double          mat_sprod(Vec * v, Mat * a, Vec * w)
212 {
213         int             i, j;
214         double          sum = 0;
215         if (a == NULL || v == NULL)
216                 return (0);
217         for (i = 0; i < a->m; i++)
218                 for (j = 0; j < a->n; j++)
219                         sum += v->el[i] * a->el[i][j] * w->el[j];
220         return (sum);
221 }
222 
223 void            mat_sum(Mat * a, Mat * b)
224 {
225         int             i, j;
226         if (a == NULL || b == NULL)
227                 return;
228         for (i = 0; i < a->m; i++)
229                 for (j = 0; j < a->n; j++)
230                         a->el[i][j] += b->el[i][j];
231 }
232 
233 void            mat_diff(Mat * a, Mat * b)
234 {
235         int             i, j;
236         if (a == NULL || b == NULL)
237                 return;
238         for (i = 0; i < a->m; i++)
239                 for (j = 0; j < a->n; j++)
240                         a->el[i][j] -= b->el[i][j];
241 }
242 
243 void            mat_times(double k, Mat * a)
244 {
245         int             i, j;
246         if (a == NULL)
247                 return;
248         for (i = 0; i < a->m; i++)
249                 for (j = 0; j < a->n; j++)
250                         a->el[i][j] *= k;
251 }
252 
253 void            mat_minus(Mat * a)
254 {
255         int             i, j;
256         if (a == NULL)
257                 return;
258         for (i = 0; i < a->m; i++)
259                 for (j = 0; j < a->n; j++)
260                         a->el[i][j] = -a->el[i][j];
261 }
262 
263 void            mat_accum(Mat * a, double k, Mat * b)
264 {
265         int             i, j;
266         if (a == NULL || b == NULL)
267                 return;
268         for (i = 0; i < a->m; i++)
269                 for (j = 0; j < a->n; j++)
270                         a->el[i][j] += k * b->el[i][j];
271 }
272 
273 void            mat_sum_tensor(Mat * a, Vec * v, Vec * w)
274 {
275         int             i, j;
276         if (a == NULL || v == NULL || w == NULL)
277                 return;
278         for (i = 0; i < a->m; i++)
279                 for (j = 0; j < a->n; j++)
280                         a->el[i][j] += v->el[i] * w->el[j];
281 }
282 
283 void            mat_accum_tensor(Mat * a, double k, Vec * v, Vec * w)
284 {
285         int             i, j;
286         if (a == NULL || v == NULL || w == NULL)
287                 return;
288         for (i = 0; i < a->m; i++)
289                 for (j = 0; j < a->n; j++)
290                         a->el[i][j] += k * v->el[i] * w->el[j];
291 }
292 
293 void            mat_print(FILE * fp, char *fmt, Mat * a)
294 {
295         int             i, j;
296         if (a == NULL)
297         {
298                 fprintf(fp, "0 0 :\n");
299                 return;
300         }
301         fprintf(fp, "%d %d :\n", a->m, a->n);
302         for (i = 0; i < a->m; i++)
303         {
304                 for (j = 0; j < a->n; j++)
305                         fprintf(fp, fmt, a->el[i][j]);
306                 fprintf(fp, "\n");
307         }
308 }
309 
310 Mat            *mat_read(FILE * fp)
311 {
312         Mat            *a;
313         int             i, j, m, n;
314         fscanf(fp, "%d %d %*s", &m, &n);
315         a = mat_make(m, n);
316         for (i = 0; i < a->m; i++)
317                 for (j = 0; j < a->n; j++)
318                         fscanf(fp, "%lf", &a->el[i][j]);
319         return (a);
320 }
321 
322 void            mat_row(Vec * row, Mat * a, int i)
323 {
324         int             j;
325         if (row == NULL || a == NULL)
326                 return;
327         for (j = 0; j < a->n; j++)
328                 row->el[j] = a->el[i][j];
329 }
330 
331 void            mat_row_get(Mat * a, int i, Vec * row)
332 {
333         int             j;
334         if (row == NULL || a == NULL)
335                 return;
336         for (j = 0; j < row->n; j++)
337                 row->el[j] = a->el[i][j];
338 }
339 
340 void            mat_row_set(Mat * a, int i, Vec * row)
341 {
342         int             j;
343         if (row == NULL || a == NULL)
344                 return;
345         for (j = 0; j < row->n; j++)
346                 a->el[i][j] = row->el[j];
347 }
348 
349 void            mat_col_get(Mat * a, int j, Vec * col)
350 {
351         int             i;
352         if (col == NULL || a == NULL)
353                 return;
354         for (i = 0; i < col->n; i++)
355                 col->el[i] = a->el[i][j];
356 }
357 
358 void            mat_col_set(Mat * a, int j, Vec * col)
359 {
360         int             i;
361         if (col == NULL || a == NULL)
362                 return;
363         for (i = 0; i < col->n; i++)
364                 a->el[i][j] = col->el[i];
365 }
366 
367 void            mat_block_get(Mat * a, int li, int lj, Mat * b)
368 {
369         int             i, j;
370         for (i = 0; i < b->m; i++)
371                 for (j = 0; j < b->n; j++)
372                         b->el[i][j] = a->el[i + li][j + lj];
373 }
374 
375 void            mat_block_set(Mat * a, int li, int lj, Mat * b)
376 {
377         int             i, j;
378         for (i = 0; i < b->m; i++)
379                 for (j = 0; j < b->n; j++)
380                         a->el[i + li][j + lj] = b->el[i][j];
381 }
382 
383 void            mat_index_get(Mat * a, Ivec * indi, Ivec * indj, Mat * b)
384 {
385         int             i, j;
386         for (i = 0; i < b->m; i++)
387                 for (j = 0; j < b->n; j++)
388                         b->el[i][j] = a->el[indi->el[i]][indj->el[j]];
389 }
390 
391 void            mat_index_set(Mat * a, Ivec * indi, Ivec * indj, Mat * b)
392 {
393         int             i, j;
394         for (i = 0; i < b->m; i++)
395                 for (j = 0; j < b->n; j++)
396                         a->el[indi->el[i]][indj->el[j]] = b->el[i][j];
397 }
398 
399 void            mat_mat3(Mat * a, Mat3 a3)
400 {
401         a->el[0][0] = mat3_xx(a3);
402         a->el[0][1] = mat3_xy(a3);
403         a->el[0][2] = mat3_xz(a3);
404         a->el[1][0] = mat3_yx(a3);
405         a->el[1][1] = mat3_yy(a3);
406         a->el[1][2] = mat3_yz(a3);
407         a->el[2][0] = mat3_zx(a3);
408         a->el[2][1] = mat3_zy(a3);
409         a->el[2][2] = mat3_zz(a3);
410 }
411 
412 Mat3            mat3_mat(Mat * a)
413 {
414         Mat3            a3;
415         mat3_xx(a3) = a->el[0][0];
416         mat3_xy(a3) = a->el[0][1];
417         mat3_xz(a3) = a->el[0][2];
418         mat3_yx(a3) = a->el[1][0];
419         mat3_yy(a3) = a->el[1][1];
420         mat3_yz(a3) = a->el[1][2];
421         mat3_zx(a3) = a->el[2][0];
422         mat3_zy(a3) = a->el[2][1];
423         mat3_zz(a3) = a->el[2][2];
424         a3.ts_id = Mat3_id;
425         return (a3);
426 }
427 

~ [ source navigation ] ~ [ diff markup ] ~ [ identifier search ] ~ [ freetext search ] ~ [ file search ] ~

This page was automatically generated by the LXR engine.
Visit the LXR main site for more information.