PageRenderTime 63ms CodeModel.GetById 38ms RepoModel.GetById 0ms app.codeStats 0ms

/SiemensSuite/TotInfo/TotInfo/Version3.cs

https://github.com/lukesandberg/PexFaultLocalization
C# | 381 lines | 276 code | 72 blank | 33 comment | 40 complexity | e55d805bd6e0d511a97cb8a1a4df8cbf MD5 | raw file
  1. using System;
  2. using System.Collections.Generic;
  3. using System.Linq;
  4. using System.Text;
  5. using System.IO;
  6. using System.Text.RegularExpressions;
  7. using Edu.Nlu.Sir.Siemens.Shared;
  8. namespace Edu.Nlu.Sir.Siemens.TotInfo
  9. {
  10. public class BaseVersion : ITotInfo, FaultyVersion
  11. {
  12. public int[] FaultLines { get { return new int[] { -1 }; } }
  13. public FaultType FaultType { get { return FaultType.NO_FAULT; } }
  14. public const int MAXLINE = 256;
  15. public const int MAXTBL = 1000;
  16. public const char COMMENT = '#';
  17. public const int NULL = 0;
  18. public const int EXIT_FAILURE = -1;
  19. public const int EXIT_SUCCESS = 0;
  20. public const int ITMAX = 100;
  21. public const double EPS = 3.0e-7;
  22. /// <summary>
  23. /// row/column header input line */
  24. /// </summary>
  25. public char[] line = new char[MAXLINE];
  26. /// <summary>
  27. /// frequency tallies
  28. /// </summary>
  29. public long[] f = new long[MAXTBL];
  30. /// <summary>
  31. /// # of rows
  32. /// </summary>
  33. public int r;
  34. /// <summary>
  35. /// # of columns
  36. /// </summary>
  37. public int c;
  38. void Main(string[] args)
  39. {
  40. char p; /* input line scan location */
  41. int i; /* row index */
  42. int j; /* column index */
  43. double info; /* computed information measure */
  44. int infodf; /* degrees of freedom for information */
  45. double totinfo = 0.0; /* accumulated information */
  46. int totdf; /* accumulated degrees of freedom */
  47. totdf = 0;
  48. string lineStr;
  49. while ((lineStr = Console.ReadLine()) != null) /* start new table */
  50. {
  51. line = lineStr.ToCharArray();
  52. p = line[0];
  53. for (int pi = 1; p != '\0' && Char.IsWhiteSpace(p); ++pi)
  54. {
  55. p = line[pi];
  56. }
  57. if (p == '\0')
  58. continue; /* skip blank line */
  59. if (p == COMMENT)
  60. { /* copy comment through */
  61. Console.WriteLine(line);
  62. continue;
  63. }
  64. string[] tokens = p.ToString().Split(" ".ToCharArray());
  65. bool parsable = Int32.TryParse(tokens[0], out r) && Int32.TryParse(tokens[1], out c);
  66. if (!parsable)
  67. {
  68. Console.WriteLine("* invalid row/column line *\n");
  69. Environment.Exit(EXIT_FAILURE);
  70. }
  71. if (r * c > MAXTBL)
  72. {
  73. Console.WriteLine("* table too large *\n");
  74. Environment.Exit(EXIT_FAILURE);
  75. }
  76. /* input tallies */
  77. string valuesString = Console.In.ReadToEnd();
  78. string[] valuesLine = valuesString.Split("\n".ToCharArray());
  79. string[][] values = (string[][])valuesLine.Select(l => l.Split(" ".ToCharArray()));
  80. long value;
  81. bool parsed;
  82. for (i = 0; i < r; ++i)
  83. {
  84. for (j = 0; j < c; ++j)
  85. {
  86. parsed = Int64.TryParse(values[i][j], out value);
  87. if (parsed)
  88. {
  89. f[i * c + j] = value;
  90. }
  91. else
  92. {
  93. Console.Write("* EOF in table *\n");
  94. Environment.Exit(EXIT_FAILURE);
  95. }
  96. }
  97. }
  98. /* compute statistic */
  99. info = InfoTbl(r, c, f,out infodf);
  100. /* print results */
  101. if (info >= 0.0)
  102. {
  103. Console.Write(String.Format("2info = {0,000.00}\tdf = {1,##}\tq = {2,000.0000}%7.4f\n",
  104. info, infodf,
  105. QChiSq(info, infodf))
  106. );
  107. totinfo += info;
  108. totdf += infodf;
  109. }
  110. else
  111. Console.Write(info < -3.5 ? "out of memory\n"
  112. : info < -2.5 ? "table too small\n"
  113. : info < -1.5 ? "negative freq\n"
  114. : "table all zeros\n");
  115. }
  116. if (totdf <= 0)
  117. {
  118. Console.Write("\n*** no information accumulated ***\n");
  119. Environment.Exit(EXIT_FAILURE);
  120. }
  121. Console.Write("\ntotal 2info = %5.2f\tdf = %2d\tq = %7.4f\n",
  122. totinfo, totdf,
  123. QChiSq(totinfo, totdf)
  124. );
  125. Environment.Exit(EXIT_SUCCESS);
  126. }
  127. public override double LGamma(double x)
  128. {
  129. double[] cof =
  130. {
  131. 76.18009173, -86.50532033, 24.01409822,
  132. -1.231739516, 0.120858003e-2, -0.536382e-5
  133. };
  134. double tmp, ser;
  135. int j;
  136. if (--x < 0.0) /* use reflection formula for accuracy */
  137. {
  138. double pix = Math.PI * x;
  139. return Math.Log(pix / Math.Sin(pix)) - LGamma(1.0 - x);
  140. }
  141. tmp = x + 5.5;
  142. tmp -= (x + 0.5) * Math.Log(tmp);
  143. ser = 1.0;
  144. for (j = 0; j < 6; ++j)
  145. ser += cof[j] / ++x;
  146. return -tmp + Math.Log(2.50662827465 * ser);
  147. }
  148. public override double
  149. gser(double a, double x)
  150. {
  151. double ap, del, sum;
  152. int n;
  153. if (x <= 0.0)
  154. return 0.0;
  155. del = sum = 1.0 / (ap = a);
  156. for (n = 1; n <= ITMAX; ++n)
  157. {
  158. sum += del *= x / ++ap;
  159. if (Math.Abs(del) < Math.Abs(sum) * EPS)
  160. return sum * Math.Exp(-x + a * Math.Log(x) - LGamma(a));
  161. }
  162. throw new ApplicationException("Execution should not have reached this line");
  163. }
  164. public override double
  165. gcf(double a, double x)
  166. {
  167. int n;
  168. double gold = 0.0, fac = 1.0, b1 = 1.0,
  169. b0 = 0.0, a0 = 1.0, a1 = x;
  170. for (n = 1; n <= ITMAX; ++n)
  171. {
  172. double anf;
  173. double an = (double)n;
  174. double ana = an - a;
  175. a0 = (a1 + a0 * ana) * fac;
  176. b0 = (b1 + b0 * ana) * fac;
  177. anf = an * fac;
  178. b1 = x * b0 + anf * b1;
  179. a1 = x * a0 + anf * a1;
  180. if (a1 != 0.0)
  181. { /* renormalize */
  182. double g = b1 * (fac = 1.0 / a1);
  183. gold = g - gold;
  184. if (Math.Abs(gold) < EPS * Math.Abs(g))
  185. return Math.Exp(-x + a * Math.Log(x) - LGamma(a)) * g;
  186. gold = g;
  187. }
  188. }
  189. throw new ApplicationException("Execution should not have reached this line");
  190. }
  191. public override double
  192. QGamma(double a, double x)
  193. {
  194. return x < a + 1.0 ? 1.0 - gser(a, x) : gcf(a, x);
  195. }
  196. public override double
  197. QChiSq(double chisq, int df)
  198. {
  199. return QGamma((double)df / 2.0, chisq / 2.0);
  200. }
  201. /// <summary>
  202. /// InfoTbl -- Kullback's information measure for a 2-way contingency table
  203. ///last edit: 88/09/19 D A Gwyn
  204. ///SCCS ID: @(#)info.c 1.1 (edited for publication)
  205. ///Special return values:
  206. ///-1.0 entire table consisted of 0 entries
  207. ///-2.0 invalid table entry (frequency less than 0)
  208. ///-3.0 invalid table dimensions (r or c less than 2)
  209. ///-4.0 unable to allocate enough working storage
  210. /// </summary>
  211. /// <param name="r"># rows in table</param>
  212. /// <param name="c"># columns in table </param>
  213. /// <param name="f">r*c frequency tallies</param>
  214. /// <param name="pdf">return # d.f. for chi-square</param>
  215. /// <returns></returns>
  216. public override double
  217. InfoTbl(int r, int c, long[] f, out int pdf)
  218. {
  219. int i; /* row index */
  220. int j; /* column index */
  221. double N; /* (double)n */
  222. double info; /* accumulates information measure */
  223. double[] xi; /* row sums */
  224. double[] xj; /* col sums */
  225. int rdf = r - 1; /* row degrees of freedom */
  226. int cdf = c - 1; /* column degrees of freedom */
  227. if (rdf <= 0 || cdf <= 0)
  228. {
  229. pdf = 0;
  230. info = -3.0;
  231. goto ret3;
  232. }
  233. pdf = rdf * cdf; /* total degrees of freedom */
  234. try {
  235. xi = new double[r];
  236. }
  237. catch(OutOfMemoryException)
  238. {
  239. info = -4.0;
  240. goto ret3;
  241. }
  242. try {
  243. xj = new double[c];
  244. }
  245. catch(OutOfMemoryException)
  246. {
  247. info = -4.0;
  248. goto ret2;
  249. }
  250. /* compute row sums and total */
  251. N = 0.0;
  252. for (i = 0; i < r; ++i)
  253. {
  254. double sum = 0.0; /* accumulator */
  255. for (j = 0; j < c; ++j)
  256. {
  257. long k = f[(i) * c + (j)];
  258. if (k < 0L)
  259. {
  260. info = -2.0;
  261. goto ret1;
  262. }
  263. sum += (double)k;
  264. }
  265. N += xi[i] = sum;
  266. }
  267. if (N <= 0.0)
  268. {
  269. info = -1.0;
  270. goto ret1;
  271. }
  272. /* compute column sums */
  273. for (j = 0; j < c; ++j)
  274. {
  275. double sum = 0.0; /* accumulator */
  276. for (i = 0; i < r; ++i)
  277. sum += (double)f[(i) * c + (j)];
  278. xj[j] = sum;
  279. }
  280. /* compute information measure (four parts) */
  281. info = N * Math.Log(N); /* part 1 */
  282. for (i = 0; i < r; ++i)
  283. {
  284. double pi = xi[i]; /* row sum */
  285. if (pi > 0.0)
  286. info -= pi * Math.Log(pi); /* part 2 */
  287. for (j = 0; j < c; ++j)
  288. {
  289. double pij = (double)f[(i) * c + (j)];
  290. if (pij > 0.0)
  291. info += pij * Math.Log(pij); /* part 3 */
  292. }
  293. }
  294. for (j = 0; j < c; ++j)
  295. {
  296. double pj = xj[j]; /* column sum */
  297. if (pj > 0.0)
  298. info -= pj * Math.Log(pj); /* part 4 */
  299. }
  300. info *= 2.0; /* for comparability with chi-square */
  301. ret1:
  302. xj = null;
  303. ret2:
  304. xi = null;
  305. ret3:
  306. return info;
  307. }
  308. }
  309. }