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 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 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 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 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 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 88 inline int db_MaxIndex2(double s[2]) 89 { 90 if(s[0]>=s[1]) return(0); 91 return(1); 92 } 93 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 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 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 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 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 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 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 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 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 */ 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