1 /****************************************************************************/
2 /*                                 matrix.c                                 */
3 /****************************************************************************/
4 /*                                                                          */
5 /* type MATRIX                                                              */
6 /*                                                                          */
7 /* Copyright (C) 1992-1995 Tomas Skalicky. All rights reserved.             */
8 /*                                                                          */
9 /****************************************************************************/
10 /*                                                                          */
11 /*        ANY USE OF THIS CODE CONSTITUTES ACCEPTANCE OF THE TERMS          */
12 /*              OF THE COPYRIGHT NOTICE (SEE FILE COPYRGHT.H)               */
13 /*                                                                          */
14 /****************************************************************************/
15 
16 #include <stddef.h>
17 #include <stdlib.h>
18 #include <math.h>
19 #include <string.h>
20 
21 #include "laspack/matrix.h"
22 #include "laspack/errhandl.h"
23 #include "laspack/copyrght.h"
24 
25 static ElType ZeroEl = { 0, 0.0 };
26 
27 static int ElCompar(const void *El1, const void *El2);
28 
M_Constr(Matrix * M,char * Name,size_t RowDim,size_t ClmDim,ElOrderType ElOrder,InstanceType Instance,Boolean OwnData)29 void M_Constr(Matrix *M, char *Name, size_t RowDim, size_t ClmDim,
30               ElOrderType ElOrder, InstanceType Instance, Boolean OwnData)
31 /* constructor of the type Matrix */
32 {
33     size_t Dim, RoC;
34 
35     M->Name = (char *)malloc((strlen(Name) + 1) * sizeof(char));
36     if (M->Name != NULL)
37         strcpy(M->Name, Name);
38     else
39         LASError(LASMemAllocErr, "M_Constr", Name, NULL, NULL);
40     M->RowDim = RowDim;
41     M->ClmDim = ClmDim;
42     M->ElOrder = ElOrder;
43     M->Instance = Instance;
44     M->LockLevel = 0;
45     M->Multipl = 1.0;
46     M->OwnData = OwnData;
47     if (OwnData) {
48         if (LASResult() == LASOK) {
49             if (ElOrder == Rowws)
50                 Dim = RowDim;
51             else
52                 Dim = ClmDim;
53 	    M->Len = (size_t *)malloc((Dim + 1) * sizeof(size_t));
54 	    M->El = (ElType **)malloc((Dim + 1) * sizeof(ElType *));
55 	    M->ElSorted = (Boolean *)malloc(sizeof(Boolean));
56 	    if (M->Len != NULL && M->El != NULL) {
57                 for (RoC = 1; RoC <= Dim; RoC++) {
58                     M->Len[RoC] = 0;
59                     M->El[RoC] = NULL;
60                 }
61                 *M->ElSorted = False;
62             } else {
63 	        LASError(LASMemAllocErr, "M_Constr", Name, NULL, NULL);
64             }
65         } else {
66 	    M->Len = NULL;
67 	    M->El = NULL;
68 	    M->ElSorted = NULL;
69         }
70     }
71 }
72 
M_Destr(Matrix * M)73 void M_Destr(Matrix *M)
74 /* destructor of the type Matrix */
75 {
76     size_t Dim, RoC;
77 
78     if (M->Name != NULL)
79         free(M->Name);
80     if (M->ElOrder == Rowws)
81         Dim = M->RowDim;
82     else
83         Dim = M->ClmDim;
84     if (M->OwnData) {
85 	if (M->Len != NULL && M->El != NULL) {
86             for (RoC = 1; RoC <= Dim; RoC++) {
87                 if (M->Len[RoC] > 0) {
88                     if (M->El[RoC] != NULL)
89                         free(M->El[RoC]);
90                 }
91             }
92         }
93         if (M->Len != NULL) {
94             free(M->Len);
95             M->Len = NULL;
96         }
97         if (M->El != NULL) {
98             free(M->El);
99             M->El = NULL;
100         }
101         if (M->ElSorted != NULL) {
102             free(M->ElSorted);
103             M->ElSorted = NULL;
104         }
105     }
106 }
107 
M_SetName(Matrix * M,char * Name)108 void M_SetName(Matrix *M, char *Name)
109 /* (re)set name of the matrix M */
110 {
111     if (LASResult() == LASOK) {
112         free(M->Name);
113         M->Name = (char *)malloc((strlen(Name) + 1) * sizeof(char));
114         if (M->Name != NULL)
115             strcpy(M->Name, Name);
116         else
117             LASError(LASMemAllocErr, "M_SetName", Name, NULL, NULL);
118     }
119 }
120 
M_GetName(Matrix * M)121 char *M_GetName(Matrix *M)
122 /* returns the name of the matrix M */
123 {
124     if (LASResult() == LASOK)
125         return(M->Name);
126     else
127         return("");
128 }
129 
M_GetRowDim(Matrix * M)130 size_t M_GetRowDim(Matrix *M)
131 /* returns the row dimension of the matrix M */
132 {
133     size_t Dim;
134 
135     if (LASResult() == LASOK)
136         Dim = M->RowDim;
137     else
138         Dim = 0;
139     return(Dim);
140 }
141 
M_GetClmDim(Matrix * M)142 size_t M_GetClmDim(Matrix *M)
143 /* returns the column dimension of the matrix M */
144 {
145     size_t Dim;
146 
147     if (LASResult() == LASOK)
148         Dim = M->ClmDim;
149     else
150         Dim = 0;
151     return(Dim);
152 }
153 
M_GetElOrder(Matrix * M)154 ElOrderType M_GetElOrder(Matrix *M)
155 /* returns the element order */
156 {
157     ElOrderType ElOrder;
158 
159     if (LASResult() == LASOK) {
160         ElOrder = M->ElOrder;
161     } else {
162         ElOrder = (ElOrderType)0;
163     }
164     return(ElOrder);
165 }
166 
M_SetLen(Matrix * M,size_t RoC,size_t Len)167 void M_SetLen(Matrix *M, size_t RoC, size_t Len)
168 /* set the length of a row or column of the matrix M */
169 {
170     size_t ElCount;
171     ElType *PtrEl;
172 
173     if (LASResult() == LASOK) {
174         if (M->Instance == Normal
175             && ((M->ElOrder == Rowws && RoC > 0 && RoC <= M->RowDim)
176             || (M->ElOrder == Clmws && RoC > 0 && RoC <= M->ClmDim))) {
177             M->Len[RoC] = Len;
178 
179             PtrEl = M->El[RoC];
180 
181             if (PtrEl != NULL) {
182                 free(PtrEl);
183 		PtrEl = NULL;
184             }
185 
186             if (Len > 0) {
187                 PtrEl = (ElType *)malloc(Len * sizeof(ElType));
188                 M->El[RoC] = PtrEl;
189 
190                 if (PtrEl != NULL) {
191                     for (ElCount = Len; ElCount > 0; ElCount--) {
192                         *PtrEl = ZeroEl;
193                         PtrEl++;
194                     }
195                 } else {
196                     LASError(LASMemAllocErr, "M_SetLen", M->Name, NULL, NULL);
197                 }
198             } else {
199                 M->El[RoC] = NULL;
200             }
201         } else {
202             if (M->Instance == Normal)
203                 LASError(LASLValErr, "M_SetLen", M->Name, NULL, NULL);
204             else
205                 LASError(LASRangeErr, "M_SetLen", M->Name, NULL, NULL);
206         }
207     }
208 }
209 
M_GetLen(Matrix * M,size_t RoC)210 size_t M_GetLen(Matrix *M, size_t RoC)
211 /* returns the length of a row or column of the matrix M */
212 {
213     size_t Len;
214 
215     if (LASResult() == LASOK) {
216         if ((M->ElOrder == Rowws && RoC > 0 && RoC <= M->RowDim) ||
217             (M->ElOrder == Clmws && RoC > 0 && RoC <= M->ClmDim)) {
218             Len = M->Len[RoC];
219         } else {
220             LASError(LASRangeErr, "M_GetLen", M->Name, NULL, NULL);
221             Len = 0;
222         }
223     } else {
224         Len = 0;
225     }
226     return(Len);
227 }
228 
M_SetEntry(Matrix * M,size_t RoC,size_t Entry,size_t Pos,Real Val)229 void M_SetEntry(Matrix *M, size_t RoC, size_t Entry, size_t Pos, Real Val)
230 /* set a new matrix entry */
231 {
232     if (LASResult() == LASOK) {
233         if ((M->ElOrder == Rowws && RoC > 0 && RoC <= M->RowDim && Pos > 0 && Pos <= M->ClmDim) ||
234             ((M->ElOrder == Clmws && RoC > 0 && RoC <= M->ClmDim && Pos > 0 && Pos <= M->RowDim) &&
235             (Entry < M->Len[RoC]))) {
236             M->El[RoC][Entry].Val = Val;
237             M->El[RoC][Entry].Pos = Pos;
238         } else {
239             LASError(LASRangeErr, "M_SetEntry", M->Name, NULL, NULL);
240         }
241     }
242 }
243 
M_GetPos(Matrix * M,size_t RoC,size_t Entry)244 size_t M_GetPos(Matrix *M, size_t RoC, size_t Entry)
245 /* returns the position of a matrix entry */
246 {
247     size_t Pos;
248 
249     if (LASResult() == LASOK)
250         if ((M->ElOrder == Rowws && RoC > 0 && RoC <= M->RowDim) ||
251             ((M->ElOrder == Clmws && RoC > 0 && RoC <= M->ClmDim) &&
252             (Entry < M->Len[RoC]))) {
253             Pos = M->El[RoC][Entry].Pos;
254         } else {
255             LASError(LASRangeErr, "M_GetPos", M->Name, NULL, NULL);
256 	    Pos = 0;
257         }
258     else
259         Pos = 0;
260     return(Pos);
261 }
262 
M_GetVal(Matrix * M,size_t RoC,size_t Entry)263 Real M_GetVal(Matrix *M, size_t RoC, size_t Entry)
264 /* returns the value of a matrix entry */
265 {
266     Real Val;
267 
268     if (LASResult() == LASOK)
269         if ((M->ElOrder == Rowws && RoC > 0 && RoC <= M->RowDim) ||
270             ((M->ElOrder == Clmws && RoC > 0 && RoC <= M->ClmDim) &&
271             (Entry < M->Len[RoC]))) {
272             Val = M->El[RoC][Entry].Val;
273         } else {
274             LASError(LASRangeErr, "M_GetVal", M->Name, NULL, NULL);
275 	    Val = 0.0;
276 	}
277     else
278         Val = 0.0;
279     return(Val);
280 }
281 
M_AddVal(Matrix * M,size_t RoC,size_t Entry,Real Val)282 void M_AddVal(Matrix *M, size_t RoC, size_t Entry, Real Val)
283 /* add a value to a matrix entry */
284 {
285     if (LASResult() == LASOK) {
286         if ((M->ElOrder == Rowws && RoC > 0 && RoC <= M->RowDim) ||
287             ((M->ElOrder == Clmws && RoC > 0 && RoC <= M->ClmDim) &&
288             (Entry < M->Len[RoC])))
289             M->El[RoC][Entry].Val += Val;
290         else
291             LASError(LASRangeErr, "M_AddVal", M->Name, NULL, NULL);
292     }
293 }
294 
M_GetEl(Matrix * M,size_t Row,size_t Clm)295 Real M_GetEl(Matrix *M, size_t Row, size_t Clm)
296 /* returns the value of a matrix element (all matrix elements are considered) */
297 {
298     Real Val;
299 
300     size_t Len, ElCount;
301     ElType *PtrEl;
302 
303     if (LASResult() == LASOK) {
304         if (Row > 0 && Row <= M->RowDim && Clm > 0 && Clm <= M->ClmDim) {
305             Val = 0.0;
306             if (M->ElOrder == Rowws) {
307                 Len = M->Len[Row];
308                 PtrEl = M->El[Row];
309                 for (ElCount = Len; ElCount > 0; ElCount--) {
310                     if ((*PtrEl).Pos == Clm)
311                         Val = (*PtrEl).Val;
312                     PtrEl++;
313                 }
314             } else if (M->ElOrder == Clmws) {
315                 Len = M->Len[Clm];
316                 PtrEl = M->El[Clm];
317                 for (ElCount = Len; ElCount > 0; ElCount--) {
318                     if ((*PtrEl).Pos == Row)
319                         Val = (*PtrEl).Val;
320                     PtrEl++;
321                 }
322             }
323         } else {
324             LASError(LASRangeErr, "M_GetEl", M->Name, NULL, NULL);
325             Val = 0.0;
326         }
327     } else {
328         Val = 0.0;
329     }
330     return(Val);
331 }
332 
M_SortEl(Matrix * M)333 void M_SortEl(Matrix *M)
334 /* sorts elements of a row or column in ascended order */
335 {
336     size_t Dim = 0, RoC;
337 
338     if (LASResult() == LASOK && !(*M->ElSorted)) {
339         if (M->ElOrder == Rowws)
340 	    Dim = M->ClmDim;
341         if (M->ElOrder == Clmws)
342 	    Dim = M->ClmDim;
343         for (RoC = 1; RoC <= Dim; RoC++) {
344             /* sort of elements by the quick sort algorithms */
345             qsort((void *)M->El[RoC], M->Len[RoC], sizeof(ElType), ElCompar);
346         }
347 
348         *M->ElSorted = True;
349     }
350 }
351 
ElCompar(const void * El1,const void * El2)352 static int ElCompar(const void *El1, const void *El2)
353 /* compares positions of two matrix elements */
354 {
355     int Compar;
356 
357     Compar = 0;
358     if (((ElType *)El1)->Pos < ((ElType *)El2)->Pos)
359         Compar = -1;
360     if (((ElType *)El1)->Pos > ((ElType *)El2)->Pos)
361         Compar = +1;
362 
363     return(Compar);
364 }
365 
M_Lock(Matrix * M)366 void M_Lock(Matrix *M)
367 /* lock the matrix M */
368 {
369     if (M != NULL)
370         M->LockLevel++;
371 }
372 
M_Unlock(Matrix * M)373 void M_Unlock(Matrix *M)
374 /* unlock the matrix M */
375 {
376     if (M != NULL) {
377         M->LockLevel--;
378         if (M->Instance == Tempor && M->LockLevel <= 0) {
379             M_Destr(M);
380 	    free(M);
381 	}
382     }
383 }
384