1 /*
2 * Copyright (C) 2011 The Android Open Source Project
3 *
4 * Licensed under the Apache License, Version 2.0 (the "License");
5 * you may not use this file except in compliance with the License.
6 * You may obtain a copy of the License at
7 *
8 * http://www.apache.org/licenses/LICENSE-2.0
9 *
10 * Unless required by applicable law or agreed to in writing, software
11 * distributed under the License is distributed on an "AS IS" BASIS,
12 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13 * See the License for the specific language governing permissions and
14 * limitations under the License.
15 */
16
17 /* $Id: db_utilities_indexing.h,v 1.3 2011/06/17 14:03:31 mbansal Exp $ */
18
19 #ifndef DB_UTILITIES_INDEXING
20 #define DB_UTILITIES_INDEXING
21
22
23
24 /*****************************************************************
25 * Lean and mean begins here *
26 *****************************************************************/
27
28 #include "db_utilities.h"
29
30 /*!
31 * \defgroup LMIndexing (LM) Indexing Utilities (Order Statistics, Matrix Operations)
32 */
33 /*\{*/
34
db_SetupMatrixRefs(double ** ar,long rows,long cols,double * a)35 inline void db_SetupMatrixRefs(double **ar,long rows,long cols,double *a)
36 {
37 long i;
38 for(i=0;i<rows;i++) ar[i]=&a[i*cols];
39 }
40
db_SymmetricExtendUpperToLower(double ** A,int rows,int cols)41 inline void db_SymmetricExtendUpperToLower(double **A,int rows,int cols)
42 {
43 int i,j;
44 for(i=1;i<rows;i++) for(j=0;j<i;j++) A[i][j]=A[j][i];
45 }
46
db_MultiplyMatrixVectorAtb(double * c,const double * const * At,const double * b,int arows,int acols)47 void inline db_MultiplyMatrixVectorAtb(double *c,const double * const *At,const double *b,int arows,int acols)
48 {
49 int i,j;
50 double acc;
51
52 for(i=0;i<arows;i++)
53 {
54 acc=0;
55 for(j=0;j<acols;j++) acc+=At[j][i]*b[j];
56 c[i]=acc;
57 }
58 }
59
db_MultiplyMatricesAB(double ** C,const double * const * A,const double * const * B,int arows,int acols,int bcols)60 inline void db_MultiplyMatricesAB(double **C,const double * const *A,const double * const *B,int arows,int acols,int bcols)
61 {
62 int i,j,k;
63 double acc;
64
65 for(i=0;i<arows;i++) for(j=0;j<bcols;j++)
66 {
67 acc=0;
68 for(k=0;k<acols;k++) acc+=A[i][k]*B[k][j];
69 C[i][j]=acc;
70 }
71 }
72
db_UpperMultiplyMatricesAtB(double ** Cu,const double * const * At,const double * const * B,int arows,int acols,int bcols)73 inline void db_UpperMultiplyMatricesAtB(double **Cu,const double * const *At,const double * const *B,int arows,int acols,int bcols)
74 {
75 int i,j,k;
76 double acc;
77
78 for(i=0;i<arows;i++) for(j=i;j<bcols;j++)
79 {
80 acc=0;
81 for(k=0;k<acols;k++) acc+=At[k][i]*B[k][j];
82 Cu[i][j]=acc;
83 }
84 }
85
86 DB_API void db_Zero(double *d,long nr);
87
db_MaxIndex2(double s[2])88 inline int db_MaxIndex2(double s[2])
89 {
90 if(s[0]>=s[1]) return(0);
91 return(1);
92 }
93
db_MaxIndex3(const double s[3])94 inline int db_MaxIndex3(const double s[3])
95 {
96 double best;
97 int pos;
98
99 best=s[0];pos=0;
100 if(s[1]>best){best=s[1];pos=1;}
101 if(s[2]>best){best=s[2];pos=2;}
102 return(pos);
103 }
104
db_MaxIndex4(const double s[4])105 inline int db_MaxIndex4(const double s[4])
106 {
107 double best;
108 int pos;
109
110 best=s[0];pos=0;
111 if(s[1]>best){best=s[1];pos=1;}
112 if(s[2]>best){best=s[2];pos=2;}
113 if(s[3]>best){best=s[3];pos=3;}
114 return(pos);
115 }
116
db_MaxIndex5(const double s[5])117 inline int db_MaxIndex5(const double s[5])
118 {
119 double best;
120 int pos;
121
122 best=s[0];pos=0;
123 if(s[1]>best){best=s[1];pos=1;}
124 if(s[2]>best){best=s[2];pos=2;}
125 if(s[3]>best){best=s[3];pos=3;}
126 if(s[4]>best){best=s[4];pos=4;}
127 return(pos);
128 }
129
db_MaxIndex6(const double s[6])130 inline int db_MaxIndex6(const double s[6])
131 {
132 double best;
133 int pos;
134
135 best=s[0];pos=0;
136 if(s[1]>best){best=s[1];pos=1;}
137 if(s[2]>best){best=s[2];pos=2;}
138 if(s[3]>best){best=s[3];pos=3;}
139 if(s[4]>best){best=s[4];pos=4;}
140 if(s[5]>best){best=s[5];pos=5;}
141 return(pos);
142 }
143
db_MaxIndex7(const double s[7])144 inline int db_MaxIndex7(const double s[7])
145 {
146 double best;
147 int pos;
148
149 best=s[0];pos=0;
150 if(s[1]>best){best=s[1];pos=1;}
151 if(s[2]>best){best=s[2];pos=2;}
152 if(s[3]>best){best=s[3];pos=3;}
153 if(s[4]>best){best=s[4];pos=4;}
154 if(s[5]>best){best=s[5];pos=5;}
155 if(s[6]>best){best=s[6];pos=6;}
156 return(pos);
157 }
158
db_MinIndex7(const double s[7])159 inline int db_MinIndex7(const double s[7])
160 {
161 double best;
162 int pos;
163
164 best=s[0];pos=0;
165 if(s[1]<best){best=s[1];pos=1;}
166 if(s[2]<best){best=s[2];pos=2;}
167 if(s[3]<best){best=s[3];pos=3;}
168 if(s[4]<best){best=s[4];pos=4;}
169 if(s[5]<best){best=s[5];pos=5;}
170 if(s[6]<best){best=s[6];pos=6;}
171 return(pos);
172 }
173
db_MinIndex9(const double s[9])174 inline int db_MinIndex9(const double s[9])
175 {
176 double best;
177 int pos;
178
179 best=s[0];pos=0;
180 if(s[1]<best){best=s[1];pos=1;}
181 if(s[2]<best){best=s[2];pos=2;}
182 if(s[3]<best){best=s[3];pos=3;}
183 if(s[4]<best){best=s[4];pos=4;}
184 if(s[5]<best){best=s[5];pos=5;}
185 if(s[6]<best){best=s[6];pos=6;}
186 if(s[7]<best){best=s[7];pos=7;}
187 if(s[8]<best){best=s[8];pos=8;}
188 return(pos);
189 }
190
db_MaxAbsIndex3(const double * s)191 inline int db_MaxAbsIndex3(const double *s)
192 {
193 double t,best;
194 int pos;
195
196 best=fabs(s[0]);pos=0;
197 t=fabs(s[1]);if(t>best){best=t;pos=1;}
198 t=fabs(s[2]);if(t>best){pos=2;}
199 return(pos);
200 }
201
db_MaxAbsIndex9(const double * s)202 inline int db_MaxAbsIndex9(const double *s)
203 {
204 double t,best;
205 int pos;
206
207 best=fabs(s[0]);pos=0;
208 t=fabs(s[1]);if(t>best){best=t;pos=1;}
209 t=fabs(s[2]);if(t>best){best=t;pos=2;}
210 t=fabs(s[3]);if(t>best){best=t;pos=3;}
211 t=fabs(s[4]);if(t>best){best=t;pos=4;}
212 t=fabs(s[5]);if(t>best){best=t;pos=5;}
213 t=fabs(s[6]);if(t>best){best=t;pos=6;}
214 t=fabs(s[7]);if(t>best){best=t;pos=7;}
215 t=fabs(s[8]);if(t>best){best=t;pos=8;}
216 return(pos);
217 }
218
219
220 /*!
221 Select ordinal pos (zero based) out of nr_elements in s.
222 temp should point to alloced memory of at least nr_elements*2
223 Optimized runtimes on 450MHz:
224 \code
225 30 with 3 microsecs
226 100 with 11 microsecs
227 300 with 30 microsecs
228 500 with 40 microsecs
229 1000 with 100 microsecs
230 5000 with 540 microsecs
231 \endcode
232 so the expected runtime is around
233 (nr_elements/10) microseconds
234 The total quickselect cost of splitting 500 hypotheses recursively
235 is thus around 100 microseconds
236
237 Does the same operation as std::nth_element().
238 */
239 DB_API double db_LeanQuickSelect(const double *s,long nr_elements,long pos,double *temp);
240
241 /*!
242 Median of 3 doubles
243 */
db_TripleMedian(double a,double b,double c)244 inline double db_TripleMedian(double a,double b,double c)
245 {
246 if(a>b)
247 {
248 if(c>a) return(a);
249 else if(c>b) return(c);
250 else return(b);
251 }
252 else
253 {
254 if(c>b) return(b);
255 else if(c>a) return(c);
256 else return(a);
257 }
258 }
259
260 /*!
261 Align float pointer to nr_bytes by moving forward
262 */
263 DB_API float* db_AlignPointer_f(float *p,unsigned long nr_bytes);
264
265 /*!
266 Align short pointer to nr_bytes by moving forward
267 */
268 DB_API short* db_AlignPointer_s(short *p,unsigned long nr_bytes);
269
270 #endif /* DB_UTILITIES_INDEXING */
271