website/archive/binaries/mac/src/muscle/intmath.h

   1 // IntMath.h: Header for doing fractional math with integers for speed.\r
   2 \r
   3 #ifndef IntMath_h\r
   4 #define IntMath_h\r
   5 \r
   6 typedef float BASETYPE;\r
   7 //typedef double BASETYPE;\r
   8 \r
   9 // Scaling factor used to store certain floating point\r
  10 // values as integers to a few significant figures.\r
  11 //const int INTSCALE = 1000;\r
  12 const int INTSCALE = 1;\r
  13 \r
  14 // Type for a probability in range 0.0 to 1.0.\r
  15 typedef BASETYPE PROB;\r
  16 \r
  17 // Type for an log-odds integer score.\r
  18 // Stored as log2(PROB)*INTSCALE.\r
  19 //typedef int   SCORE;\r
  20 typedef BASETYPE SCORE;\r
  21 \r
  22 // Type for a weight.\r
  23 // Stored as w*INTSCALE where w is in range 0.0 to 1.0.\r
  24 //typedef unsigned WEIGHT;\r
  25 typedef BASETYPE WEIGHT;\r
  26 \r
  27 // Type for a fractional weighted count stored as n*WEIGHT/N\r
  28 // where n=measured count (integer >= 0) and N is total for\r
  29 // the distribution (e.g., n=number of residues of a given\r
  30 // type in a column, N=number of residues in the column).\r
  31 // Hence values in an FCOUNT variable range from 0..INTSCALE\r
  32 // as an integer, representing "true" values 0.0 to 1.0.\r
  33 //typedef unsigned FCOUNT;\r
  34 typedef BASETYPE FCOUNT;\r
  35 \r
  36 // Representation of -infinity. Value should\r
  37 // be large and negative, but not so large\r
  38 // that adding a few of them overflows.\r
  39 // TODO: Multiplied by 10 to work around bug\r
  40 // when aligning Bali 1ckaA in ref4, which is\r
  41 // so long that B->Mmax got to -infinity, causing\r
  42 // traceback to fail.\r
  43 //const int MINUS_INFINITY = -10000000;\r
  44 const BASETYPE MINUS_INFINITY = (BASETYPE) -1e37;\r
  45 const BASETYPE PLUS_INFINITY = (BASETYPE) 1e37;\r
  46 \r
  47 // Probability relative to a null model\r
  48 typedef double RPROB;\r
  49 \r
  50 PROB ScoreToProb(SCORE Score);\r
  51 SCORE ProbToScore(PROB Prob);\r
  52 SCORE DoubleToScore(double d);\r
  53 WEIGHT DoubleToWeight(double d);\r
  54 double WeightToDouble(WEIGHT w);\r
  55 SCORE MulScoreWeight(SCORE Score, WEIGHT Weight);\r
  56 bool ScoreEq(SCORE s1, SCORE s2);\r
  57 bool BTEq(double b1, double b2);\r
  58 \r
  59 static double ScoreToDouble(SCORE Score)\r
  60         {\r
  61         return (double) Score / (double) INTSCALE;\r
  62         }\r
  63 \r
  64 #if     0\r
  65 // In-line assembler for Result = (x*y)/z\r
  66 // Note that imul and idiv will do 64-bit arithmetic\r
  67 // on 32-bit operands, so this shouldn't overflow\r
  68 // Can't write this efficiently in C/C++ (would\r
  69 // often overlow 32 bits).\r
  70 #define MulDivAssign(Result, x, y, z)   \\r
  71         {                                                                       \\r
  72         int X = (x);                                            \\r
  73         int Y = (y);                                            \\r
  74         int Z = (z);                                            \\r
  75         _asm mov        eax,X                                   \\r
  76         _asm imul       Y                                               \\r
  77         _asm mov        ecx,Z                                   \\r
  78         _asm idiv       ecx                                             \\r
  79         _asm mov        Result,eax                              \\r
  80         }\r
  81 #else\r
  82 #define MulDivAssign(Result, x, y, z)   Result = (((x)*(y))/(z))\r
  83 #endif\r
  84 \r
  85 #define MulScoreWeight(r, s, w)         MulDivAssign(r, s, w, INTSCALE)\r
  86 #define MulWeightWCount(r, wt, wc)      MulDivAssign(r, wt, wc, INTSCALE)\r
  87 #define MulFCountScore(r, fc, sc)       MulDivAssign(r, fc, sc, INTSCALE)\r
  88 \r
  89 #if     _DEBUG\r
  90 \r
  91 static inline SCORE Add2(SCORE a, SCORE b)\r
  92         {\r
  93         if (MINUS_INFINITY == a)\r
  94                 return MINUS_INFINITY;\r
  95         if (MINUS_INFINITY == b)\r
  96                 return MINUS_INFINITY;\r
  97         SCORE sum = a + b;\r
  98         if (sum < MINUS_INFINITY)\r
  99                 return MINUS_INFINITY;\r
 100 //      assert(sum < OVERFLOW_WARN);\r
 101         return sum;\r
 102         }\r
 103 \r
 104 static inline SCORE Add3(SCORE a, SCORE b, SCORE c)\r
 105         {\r
 106         return Add2(Add2(a, b), c);\r
 107         }\r
 108 \r
 109 static inline SCORE Add4(SCORE a, SCORE b, SCORE c, SCORE d)\r
 110         {\r
 111         return Add2(Add2(a, b), Add2(c, d));\r
 112         }\r
 113 \r
 114 static inline SCORE Add5(SCORE a, SCORE b, SCORE c, SCORE d, SCORE e)\r
 115         {\r
 116         return Add3(Add2(a, b), Add2(c, d), e);\r
 117         }\r
 118 \r
 119 static inline SCORE Add6(SCORE a, SCORE b, SCORE c, SCORE d, SCORE e, SCORE f)\r
 120         {\r
 121         return Add3(Add2(a, b), Add2(c, d), Add2(e, f));\r
 122         }\r
 123 \r
 124 static inline SCORE Add7(SCORE a, SCORE b, SCORE c, SCORE d, SCORE e, SCORE f, SCORE g)\r
 125         {\r
 126         return Add4(Add2(a, b), Add2(c, d), Add2(e, f), g);\r
 127         }\r
 128 \r
 129 static inline SCORE Mul2(SCORE a, SCORE b)\r
 130         {\r
 131         if (MINUS_INFINITY == a)\r
 132                 return MINUS_INFINITY;\r
 133         if (MINUS_INFINITY == b)\r
 134                 return MINUS_INFINITY;\r
 135         //__int64 prod = (__int64) a * (__int64) b;\r
 136         //assert((SCORE) prod == prod);\r
 137         //return (SCORE) prod;\r
 138         return a*b;\r
 139         }\r
 140 \r
 141 static inline SCORE Sub2(SCORE a, SCORE b)\r
 142         {\r
 143         if (MINUS_INFINITY == a)\r
 144                 return MINUS_INFINITY;\r
 145         if (MINUS_INFINITY == b)\r
 146                 return MINUS_INFINITY;\r
 147         SCORE diff = a - b;\r
 148         if (diff < MINUS_INFINITY)\r
 149                 return MINUS_INFINITY;\r
 150 //      assert(diff < OVERFLOW_WARN);\r
 151         return diff;\r
 152         }\r
 153 \r
 154 static inline SCORE Div2(SCORE a, int b)\r
 155         {\r
 156         if (MINUS_INFINITY == a)\r
 157                 return MINUS_INFINITY;\r
 158         return a/b;\r
 159         }\r
 160 \r
 161 //static inline SCORE MulScoreWeight(SCORE s, WEIGHT w)\r
 162 //      {\r
 163 //      SCORE Prod = s*(SCORE) w;\r
 164 //      assert(Prod < OVERFLOW_WARN);\r
 165 //      extern void Log(const char Format[], ...);\r
 166 //      if (Prod/(SCORE) w != s)\r
 167 //              Log("**WARRNING MulScoreWeight Prod=%d w=%d Prod/w=%d s=%d\n",\r
 168 //                Prod,\r
 169 //                w,\r
 170 //                Prod/(SCORE) w,\r
 171 //                s);\r
 172 //      assert(Prod/ (SCORE) w == s);\r
 173 //      return Prod/INTSCALE;\r
 174 //      }\r
 175 //\r
 176 //static inline WCOUNT MulWeightWCount(WEIGHT wt, WCOUNT wc)\r
 177 //      {\r
 178 //      return (wt*wc)/INTSCALE;\r
 179 //      }\r
 180 \r
 181 #else\r
 182 #define Add2(a, b)                                      ((a) + (b))\r
 183 #define Sub2(a, b)                                      ((MINUS_INFINITY == (a)) ? MINUS_INFINITY : ((a) - (b)))\r
 184 #define Div2(a, b)                                      ((MINUS_INFINITY == (a)) ? MINUS_INFINITY : ((a) / (b)))\r
 185 #define Add3(a, b, c)                           ((a) + (b) + (c))\r
 186 #define Add4(a, b, c, d)                        ((a) + (b) + (c) + (d))\r
 187 #define Add5(a, b, c, d, e)                     ((a) + (b) + (c) + (d) + (e))\r
 188 #define Add6(a, b, c, d, e, f)          ((a) + (b) + (c) + (d) + (e) + (f))\r
 189 #define Add7(a, b, c, d, e, f, g)       ((a) + (b) + (c) + (d) + (e) + (f) + (g))\r
 190 //#define       MulScoreWeight(s, w)            (((s)*(SCORE) (w))/INTSCALE)\r
 191 #define Mul2(a, b)                                      ((a)*(b))\r
 192 #endif\r
 193 \r
 194 //static inline SCORE MulFCountScore(FCOUNT fc, SCORE sc)\r
 195 //      {\r
 196 //// Fast way to say "if (fc >= 2^15 || sc >= 2^15)":\r
 197 //      if ((fc | sc) & 0xffff1000)\r
 198 //              {\r
 199 //              SCORE Score = ((fc+5)/10)*sc;\r
 200 //              assert(Score < assert);\r
 201 //              OVERFLOW_WARN(Score > MINUS_INFINITY);\r
 202 //              return Score/(INTSCALE/10);\r
 203 //              }\r
 204 //      SCORE Score = fc*sc;\r
 205 //      assert(Score < OVERFLOW_WARN);\r
 206 //      assert(Score > MINUS_INFINITY);\r
 207 //      return Score/INTSCALE;\r
 208 //      }\r
 209 \r
 210 #endif  // IntMath_h\r