1/* matrix/gsl_matrix_uchar.h
2 *
3 * Copyright (C) 1996, 1997, 1998, 1999, 2000, 2007 Gerard Jungman, Brian Gough
4 *
5 * This program is free software; you can redistribute it and/or modify
6 * it under the terms of the GNU General Public License as published by
7 * the Free Software Foundation; either version 3 of the License, or (at
8 * your option) any later version.
9 *
10 * This program is distributed in the hope that it will be useful, but
11 * WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13 * General Public License for more details.
14 *
15 * You should have received a copy of the GNU General Public License
16 * along with this program; if not, write to the Free Software
17 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
18 */
19
20#ifndef __GSL_MATRIX_UCHAR_H__
21#define __GSL_MATRIX_UCHAR_H__
22
23#include <stdlib.h>
24#include <gsl/gsl_types.h>
25#include <gsl/gsl_errno.h>
26#include <gsl/gsl_inline.h>
27#include <gsl/gsl_check_range.h>
28#include <gsl/gsl_vector_uchar.h>
29
30#undef __BEGIN_DECLS
31#undef __END_DECLS
32#ifdef __cplusplus
33# define __BEGIN_DECLS extern "C" {
34# define __END_DECLS }
35#else
36# define __BEGIN_DECLS /* empty */
37# define __END_DECLS /* empty */
38#endif
39
40__BEGIN_DECLS
41
42typedef struct
43{
44 size_t size1;
45 size_t size2;
46 size_t tda;
47 unsigned char * data;
48 gsl_block_uchar * block;
49 int owner;
50} gsl_matrix_uchar;
51
52typedef struct
53{
54 gsl_matrix_uchar matrix;
55} _gsl_matrix_uchar_view;
56
57typedef _gsl_matrix_uchar_view gsl_matrix_uchar_view;
58
59typedef struct
60{
61 gsl_matrix_uchar matrix;
62} _gsl_matrix_uchar_const_view;
63
64typedef const _gsl_matrix_uchar_const_view gsl_matrix_uchar_const_view;
65
66/* Allocation */
67
68gsl_matrix_uchar *
69gsl_matrix_uchar_alloc (const size_t n1, const size_t n2);
70
71gsl_matrix_uchar *
72gsl_matrix_uchar_calloc (const size_t n1, const size_t n2);
73
74gsl_matrix_uchar *
75gsl_matrix_uchar_alloc_from_block (gsl_block_uchar * b,
76 const size_t offset,
77 const size_t n1,
78 const size_t n2,
79 const size_t d2);
80
81gsl_matrix_uchar *
82gsl_matrix_uchar_alloc_from_matrix (gsl_matrix_uchar * m,
83 const size_t k1,
84 const size_t k2,
85 const size_t n1,
86 const size_t n2);
87
88gsl_vector_uchar *
89gsl_vector_uchar_alloc_row_from_matrix (gsl_matrix_uchar * m,
90 const size_t i);
91
92gsl_vector_uchar *
93gsl_vector_uchar_alloc_col_from_matrix (gsl_matrix_uchar * m,
94 const size_t j);
95
96void gsl_matrix_uchar_free (gsl_matrix_uchar * m);
97
98/* Views */
99
100_gsl_matrix_uchar_view
101gsl_matrix_uchar_submatrix (gsl_matrix_uchar * m,
102 const size_t i, const size_t j,
103 const size_t n1, const size_t n2);
104
105_gsl_vector_uchar_view
106gsl_matrix_uchar_row (gsl_matrix_uchar * m, const size_t i);
107
108_gsl_vector_uchar_view
109gsl_matrix_uchar_column (gsl_matrix_uchar * m, const size_t j);
110
111_gsl_vector_uchar_view
112gsl_matrix_uchar_diagonal (gsl_matrix_uchar * m);
113
114_gsl_vector_uchar_view
115gsl_matrix_uchar_subdiagonal (gsl_matrix_uchar * m, const size_t k);
116
117_gsl_vector_uchar_view
118gsl_matrix_uchar_superdiagonal (gsl_matrix_uchar * m, const size_t k);
119
120_gsl_vector_uchar_view
121gsl_matrix_uchar_subrow (gsl_matrix_uchar * m, const size_t i,
122 const size_t offset, const size_t n);
123
124_gsl_vector_uchar_view
125gsl_matrix_uchar_subcolumn (gsl_matrix_uchar * m, const size_t j,
126 const size_t offset, const size_t n);
127
128_gsl_matrix_uchar_view
129gsl_matrix_uchar_view_array (unsigned char * base,
130 const size_t n1,
131 const size_t n2);
132
133_gsl_matrix_uchar_view
134gsl_matrix_uchar_view_array_with_tda (unsigned char * base,
135 const size_t n1,
136 const size_t n2,
137 const size_t tda);
138
139
140_gsl_matrix_uchar_view
141gsl_matrix_uchar_view_vector (gsl_vector_uchar * v,
142 const size_t n1,
143 const size_t n2);
144
145_gsl_matrix_uchar_view
146gsl_matrix_uchar_view_vector_with_tda (gsl_vector_uchar * v,
147 const size_t n1,
148 const size_t n2,
149 const size_t tda);
150
151
152_gsl_matrix_uchar_const_view
153gsl_matrix_uchar_const_submatrix (const gsl_matrix_uchar * m,
154 const size_t i, const size_t j,
155 const size_t n1, const size_t n2);
156
157_gsl_vector_uchar_const_view
158gsl_matrix_uchar_const_row (const gsl_matrix_uchar * m,
159 const size_t i);
160
161_gsl_vector_uchar_const_view
162gsl_matrix_uchar_const_column (const gsl_matrix_uchar * m,
163 const size_t j);
164
165_gsl_vector_uchar_const_view
166gsl_matrix_uchar_const_diagonal (const gsl_matrix_uchar * m);
167
168_gsl_vector_uchar_const_view
169gsl_matrix_uchar_const_subdiagonal (const gsl_matrix_uchar * m,
170 const size_t k);
171
172_gsl_vector_uchar_const_view
173gsl_matrix_uchar_const_superdiagonal (const gsl_matrix_uchar * m,
174 const size_t k);
175
176_gsl_vector_uchar_const_view
177gsl_matrix_uchar_const_subrow (const gsl_matrix_uchar * m, const size_t i,
178 const size_t offset, const size_t n);
179
180_gsl_vector_uchar_const_view
181gsl_matrix_uchar_const_subcolumn (const gsl_matrix_uchar * m, const size_t j,
182 const size_t offset, const size_t n);
183
184_gsl_matrix_uchar_const_view
185gsl_matrix_uchar_const_view_array (const unsigned char * base,
186 const size_t n1,
187 const size_t n2);
188
189_gsl_matrix_uchar_const_view
190gsl_matrix_uchar_const_view_array_with_tda (const unsigned char * base,
191 const size_t n1,
192 const size_t n2,
193 const size_t tda);
194
195_gsl_matrix_uchar_const_view
196gsl_matrix_uchar_const_view_vector (const gsl_vector_uchar * v,
197 const size_t n1,
198 const size_t n2);
199
200_gsl_matrix_uchar_const_view
201gsl_matrix_uchar_const_view_vector_with_tda (const gsl_vector_uchar * v,
202 const size_t n1,
203 const size_t n2,
204 const size_t tda);
205
206/* Operations */
207
208void gsl_matrix_uchar_set_zero (gsl_matrix_uchar * m);
209void gsl_matrix_uchar_set_identity (gsl_matrix_uchar * m);
210void gsl_matrix_uchar_set_all (gsl_matrix_uchar * m, unsigned char x);
211
212int gsl_matrix_uchar_fread (FILE * stream, gsl_matrix_uchar * m) ;
213int gsl_matrix_uchar_fwrite (FILE * stream, const gsl_matrix_uchar * m) ;
214int gsl_matrix_uchar_fscanf (FILE * stream, gsl_matrix_uchar * m);
215int gsl_matrix_uchar_fprintf (FILE * stream, const gsl_matrix_uchar * m, const char * format);
216
217int gsl_matrix_uchar_memcpy(gsl_matrix_uchar * dest, const gsl_matrix_uchar * src);
218int gsl_matrix_uchar_swap(gsl_matrix_uchar * m1, gsl_matrix_uchar * m2);
219
220int gsl_matrix_uchar_swap_rows(gsl_matrix_uchar * m, const size_t i, const size_t j);
221int gsl_matrix_uchar_swap_columns(gsl_matrix_uchar * m, const size_t i, const size_t j);
222int gsl_matrix_uchar_swap_rowcol(gsl_matrix_uchar * m, const size_t i, const size_t j);
223int gsl_matrix_uchar_transpose (gsl_matrix_uchar * m);
224int gsl_matrix_uchar_transpose_memcpy (gsl_matrix_uchar * dest, const gsl_matrix_uchar * src);
225
226unsigned char gsl_matrix_uchar_max (const gsl_matrix_uchar * m);
227unsigned char gsl_matrix_uchar_min (const gsl_matrix_uchar * m);
228void gsl_matrix_uchar_minmax (const gsl_matrix_uchar * m, unsigned char * min_out, unsigned char * max_out);
229
230void gsl_matrix_uchar_max_index (const gsl_matrix_uchar * m, size_t * imax, size_t *jmax);
231void gsl_matrix_uchar_min_index (const gsl_matrix_uchar * m, size_t * imin, size_t *jmin);
232void gsl_matrix_uchar_minmax_index (const gsl_matrix_uchar * m, size_t * imin, size_t * jmin, size_t * imax, size_t * jmax);
233
234int gsl_matrix_uchar_equal (const gsl_matrix_uchar * a, const gsl_matrix_uchar * b);
235
236int gsl_matrix_uchar_isnull (const gsl_matrix_uchar * m);
237int gsl_matrix_uchar_ispos (const gsl_matrix_uchar * m);
238int gsl_matrix_uchar_isneg (const gsl_matrix_uchar * m);
239int gsl_matrix_uchar_isnonneg (const gsl_matrix_uchar * m);
240
241int gsl_matrix_uchar_add (gsl_matrix_uchar * a, const gsl_matrix_uchar * b);
242int gsl_matrix_uchar_sub (gsl_matrix_uchar * a, const gsl_matrix_uchar * b);
243int gsl_matrix_uchar_mul_elements (gsl_matrix_uchar * a, const gsl_matrix_uchar * b);
244int gsl_matrix_uchar_div_elements (gsl_matrix_uchar * a, const gsl_matrix_uchar * b);
245int gsl_matrix_uchar_scale (gsl_matrix_uchar * a, const double x);
246int gsl_matrix_uchar_add_constant (gsl_matrix_uchar * a, const double x);
247int gsl_matrix_uchar_add_diagonal (gsl_matrix_uchar * a, const double x);
248
249/***********************************************************************/
250/* The functions below are obsolete */
251/***********************************************************************/
252int gsl_matrix_uchar_get_row(gsl_vector_uchar * v, const gsl_matrix_uchar * m, const size_t i);
253int gsl_matrix_uchar_get_col(gsl_vector_uchar * v, const gsl_matrix_uchar * m, const size_t j);
254int gsl_matrix_uchar_set_row(gsl_matrix_uchar * m, const size_t i, const gsl_vector_uchar * v);
255int gsl_matrix_uchar_set_col(gsl_matrix_uchar * m, const size_t j, const gsl_vector_uchar * v);
256/***********************************************************************/
257
258/* inline functions if you are using GCC */
259
260INLINE_DECL unsigned char gsl_matrix_uchar_get(const gsl_matrix_uchar * m, const size_t i, const size_t j);
261INLINE_DECL void gsl_matrix_uchar_set(gsl_matrix_uchar * m, const size_t i, const size_t j, const unsigned char x);
262INLINE_DECL unsigned char * gsl_matrix_uchar_ptr(gsl_matrix_uchar * m, const size_t i, const size_t j);
263INLINE_DECL const unsigned char * gsl_matrix_uchar_const_ptr(const gsl_matrix_uchar * m, const size_t i, const size_t j);
264
265#ifdef HAVE_INLINE
266INLINE_FUN
267unsigned char
268gsl_matrix_uchar_get(const gsl_matrix_uchar * m, const size_t i, const size_t j)
269{
270#if GSL_RANGE_CHECK
271 if (GSL_RANGE_COND(1))
272 {
273 if (i >= m->size1)
274 {
275 GSL_ERROR_VAL("first index out of range", GSL_EINVAL, 0) ;
276 }
277 else if (j >= m->size2)
278 {
279 GSL_ERROR_VAL("second index out of range", GSL_EINVAL, 0) ;
280 }
281 }
282#endif
283 return m->data[i * m->tda + j] ;
284}
285
286INLINE_FUN
287void
288gsl_matrix_uchar_set(gsl_matrix_uchar * m, const size_t i, const size_t j, const unsigned char x)
289{
290#if GSL_RANGE_CHECK
291 if (GSL_RANGE_COND(1))
292 {
293 if (i >= m->size1)
294 {
295 GSL_ERROR_VOID("first index out of range", GSL_EINVAL) ;
296 }
297 else if (j >= m->size2)
298 {
299 GSL_ERROR_VOID("second index out of range", GSL_EINVAL) ;
300 }
301 }
302#endif
303 m->data[i * m->tda + j] = x ;
304}
305
306INLINE_FUN
307unsigned char *
308gsl_matrix_uchar_ptr(gsl_matrix_uchar * m, const size_t i, const size_t j)
309{
310#if GSL_RANGE_CHECK
311 if (GSL_RANGE_COND(1))
312 {
313 if (i >= m->size1)
314 {
315 GSL_ERROR_NULL("first index out of range", GSL_EINVAL) ;
316 }
317 else if (j >= m->size2)
318 {
319 GSL_ERROR_NULL("second index out of range", GSL_EINVAL) ;
320 }
321 }
322#endif
323 return (unsigned char *) (m->data + (i * m->tda + j)) ;
324}
325
326INLINE_FUN
327const unsigned char *
328gsl_matrix_uchar_const_ptr(const gsl_matrix_uchar * m, const size_t i, const size_t j)
329{
330#if GSL_RANGE_CHECK
331 if (GSL_RANGE_COND(1))
332 {
333 if (i >= m->size1)
334 {
335 GSL_ERROR_NULL("first index out of range", GSL_EINVAL) ;
336 }
337 else if (j >= m->size2)
338 {
339 GSL_ERROR_NULL("second index out of range", GSL_EINVAL) ;
340 }
341 }
342#endif
343 return (const unsigned char *) (m->data + (i * m->tda + j)) ;
344}
345
346#endif
347
348__END_DECLS
349
350#endif /* __GSL_MATRIX_UCHAR_H__ */
351