138 #include <config_auto.h>
142 #include "allheaders.h"
145 static const l_int32 BinSizeArray[] = {2, 5, 10, 20, 50, 100, 200, 500, 1000,\
146 2000, 5000, 10000, 20000, 50000, 100000, 200000,\
147 500000, 1000000, 2000000, 5000000, 10000000,\
148 200000000, 50000000, 100000000};
149 static const l_int32 NBinSizes = 24;
152 #ifndef NO_CONSOLE_IO
153 #define DEBUG_HISTO 0
154 #define DEBUG_CROSSINGS 0
155 #define DEBUG_FREQUENCY 0
186 l_int32 i, j, n, hsize, len;
188 l_float32 *fa, *fas, *fad;
192 return (
NUMA *)ERROR_PTR(
"nas not defined", __func__, NULL);
194 return (
NUMA *)ERROR_PTR(
"size must be > 0", __func__, NULL);
195 if ((size & 1) == 0 ) {
196 L_WARNING(
"sel size must be odd; increasing by 1\n", __func__);
210 if ((fas = (l_float32 *)LEPT_CALLOC(len,
sizeof(l_float32))) == NULL)
211 return (
NUMA *)ERROR_PTR(
"fas not made", __func__, NULL);
212 for (i = 0; i < hsize; i++)
214 for (i = hsize + n; i < len; i++)
217 for (i = 0; i < n; i++)
218 fas[hsize + i] = fa[i];
223 for (i = 0; i < n; i++) {
225 for (j = 0; j < size; j++)
226 minval = L_MIN(minval, fas[i + j]);
253 l_int32 i, j, n, hsize, len;
255 l_float32 *fa, *fas, *fad;
259 return (
NUMA *)ERROR_PTR(
"nas not defined", __func__, NULL);
261 return (
NUMA *)ERROR_PTR(
"size must be > 0", __func__, NULL);
262 if ((size & 1) == 0 ) {
263 L_WARNING(
"sel size must be odd; increasing by 1\n", __func__);
277 if ((fas = (l_float32 *)LEPT_CALLOC(len,
sizeof(l_float32))) == NULL)
278 return (
NUMA *)ERROR_PTR(
"fas not made", __func__, NULL);
279 for (i = 0; i < hsize; i++)
281 for (i = hsize + n; i < len; i++)
284 for (i = 0; i < n; i++)
285 fas[hsize + i] = fa[i];
290 for (i = 0; i < n; i++) {
292 for (j = 0; j < size; j++)
293 maxval = L_MAX(maxval, fas[i + j]);
323 return (
NUMA *)ERROR_PTR(
"nas not defined", __func__, NULL);
325 return (
NUMA *)ERROR_PTR(
"size must be > 0", __func__, NULL);
326 if ((size & 1) == 0 ) {
327 L_WARNING(
"sel size must be odd; increasing by 1\n", __func__);
365 NUMA *nab, *nat1, *nat2, *nad;
368 return (
NUMA *)ERROR_PTR(
"nas not defined", __func__, NULL);
370 return (
NUMA *)ERROR_PTR(
"size must be > 0", __func__, NULL);
371 if ((size & 1) == 0 ) {
372 L_WARNING(
"sel size must be odd; increasing by 1\n", __func__);
416 return (
NUMA *)ERROR_PTR(
"nas not defined", __func__, NULL);
419 return (
NUMA *)ERROR_PTR(
"nad not made", __func__, NULL);
421 for (i = 0; i < n; i++) {
423 val = scale * (val + shift);
450 l_float32 sum, sumsq, val, mean, var;
452 if (pmean) *pmean = 0.0;
453 if (pvar) *pvar = 0.0;
454 if (prvar) *prvar = 0.0;
455 if (!pmean && !pvar && !prvar)
456 return ERROR_INT(
"nothing requested", __func__, 1);
458 return ERROR_INT(
"na not defined", __func__, 1);
460 return ERROR_INT(
"na is empty", __func__, 1);
461 first = L_MAX(0, first);
462 if (last < 0) last = n - 1;
464 return ERROR_INT(
"invalid first", __func__, 1);
466 L_WARNING(
"last = %d is beyond max index = %d; adjusting\n",
467 __func__, last, n - 1);
471 return ERROR_INT(
"first > last\n", __func__, 1);
472 ni = last - first + 1;
474 for (i = first; i <= last; i++) {
484 var = sumsq / ni - mean * mean;
485 if (pvar) *pvar = var;
486 if (prvar) *prvar = sqrtf(var);
535 return ERROR_INT(
"nas not defined", __func__, 1);
537 L_WARNING(
"filter wider than input array!\n", __func__);
539 if (!pnav && !pnarv) {
577 l_int32 i, n, n1, width;
579 l_float32 *fa1, *fad, *suma;
583 return (
NUMA *)ERROR_PTR(
"nas not defined", __func__, NULL);
587 L_WARNING(
"filter wider than input array!\n", __func__);
596 if ((suma = (l_float32 *)LEPT_CALLOC(n1 + 1,
sizeof(l_float32))) == NULL) {
599 return (
NUMA *)ERROR_PTR(
"suma not made", __func__, NULL);
603 for (i = 0; i < n1; i++) {
608 norm = 1. / (2 * wc + 1);
609 for (i = 0; i < n; i++)
610 fad[i] = norm * (suma[width + i] - suma[i]);
635 l_int32 i, n, n1, width;
637 l_float32 *fa1, *fad, *suma;
641 return (
NUMA *)ERROR_PTR(
"nas not defined", __func__, NULL);
645 L_WARNING(
"filter wider than input array!\n", __func__);
654 if ((suma = (l_float32 *)LEPT_CALLOC(n1 + 1,
sizeof(l_float32))) == NULL) {
657 return (
NUMA *)ERROR_PTR(
"suma not made", __func__, NULL);
661 for (i = 0; i < n1; i++) {
662 sum += fa1[i] * fa1[i];
666 norm = 1. / (2 * wc + 1);
667 for (i = 0; i < n; i++)
668 fad[i] = norm * (suma[width + i] - suma[i]);
705 l_float32 *fam, *fams, *fav, *farv;
708 if (pnav) *pnav = NULL;
709 if (pnarv) *pnarv = NULL;
711 return ERROR_INT(
"neither &nav nor &narv are defined", __func__, 1);
713 return ERROR_INT(
"nam not defined", __func__, 1);
715 return ERROR_INT(
"nams not defined", __func__, 1);
719 return ERROR_INT(
"sizes of nam and nams differ", __func__, 1);
734 for (i = 0; i < nm; i++) {
735 var = fams[i] - fam[i] * fam[i];
739 farv[i] = sqrtf(var);
769 NUMA *na1, *na2, *nad;
772 return (
NUMA *)ERROR_PTR(
"nas not defined", __func__, NULL);
776 L_ERROR(
"filter too small; returning a copy\n", __func__);
780 if (halfwin > (n - 1) / 2) {
781 halfwin = (n - 1) / 2;
782 L_INFO(
"reducing filter to halfwin = %d\n", __func__, halfwin);
791 for (i = 0; i < n; i++) {
817 return (
NUMA *)ERROR_PTR(
"nas not defined", __func__, NULL);
821 return (
NUMA *)ERROR_PTR(
"nad not made", __func__, NULL);
823 for (i = 0; i < n; i++) {
866 l_int32 i, n, ival, hval;
867 l_int32 iminval, imaxval, range, binsize, nbins, ibin;
868 l_float32 val, ratio;
871 if (pbinsize) *pbinsize = 0;
872 if (pbinstart) *pbinstart = 0;
874 return (
NUMA *)ERROR_PTR(
"na not defined", __func__, NULL);
876 return (
NUMA *)ERROR_PTR(
"maxbins < 1", __func__, NULL);
880 iminval = (l_int32)(val + 0.5);
882 imaxval = (l_int32)(val + 0.5);
883 if (pbinstart == NULL) {
886 return (
NUMA *)ERROR_PTR(
"all values < 0", __func__, NULL);
890 range = imaxval - iminval + 1;
891 if (range > maxbins - 1) {
892 ratio = (l_float64)range / (l_float64)maxbins;
894 for (i = 0; i < NBinSizes; i++) {
895 if (ratio < BinSizeArray[i]) {
896 binsize = BinSizeArray[i];
901 return (
NUMA *)ERROR_PTR(
"numbers too large", __func__, NULL);
905 if (pbinsize) *pbinsize = binsize;
906 nbins = 1 + range / binsize;
909 if (pbinstart && binsize > 1) {
911 iminval = binsize * (iminval / binsize);
913 iminval = binsize * ((iminval - binsize + 1) / binsize);
915 if (pbinstart) *pbinstart = iminval;
918 lept_stderr(
" imaxval = %d, range = %d, nbins = %d\n",
919 imaxval, range, nbins);
924 return (
NUMA *)ERROR_PTR(
"nai not made", __func__, NULL);
931 return (
NUMA *)ERROR_PTR(
"nahist not made", __func__, NULL);
935 for (i = 0; i < n; i++) {
937 ibin = (ival - iminval) / binsize;
938 if (ibin >= 0 && ibin < nbins) {
975 l_int32 i, n, imin, imax, irange, ibin, ival, allints;
976 l_float32 minval, maxval, range, binsize, fval;
980 return (
NUMA *)ERROR_PTR(
"na not defined", __func__, NULL);
981 maxbins = L_MAX(1, maxbins);
992 if (allints && (maxval - minval < maxbins)) {
993 imin = (l_int32)minval;
994 imax = (l_int32)maxval;
995 irange = imax - imin + 1;
999 for (i = 0; i < n; i++) {
1010 range = maxval - minval;
1011 binsize = range / (l_float32)maxbins;
1021 for (i = 0; i < n; i++) {
1023 ibin = (l_int32)((fval - minval) / binsize);
1024 ibin = L_MIN(ibin, maxbins - 1);
1058 l_int32 i, n, nbins, ival, ibin;
1059 l_float32 val, maxval;
1063 return (
NUMA *)ERROR_PTR(
"na not defined", __func__, NULL);
1065 return (
NUMA *)ERROR_PTR(
"binsize must be > 0.0", __func__, NULL);
1066 if (binsize > maxsize)
1071 maxsize = L_MIN(maxsize, maxval);
1072 nbins = (l_int32)(maxsize / binsize) + 1;
1077 return (
NUMA *)ERROR_PTR(
"nad not made", __func__, NULL);
1080 for (i = 0; i < n; i++) {
1082 ibin = (l_int32)(val / binsize);
1083 if (ibin >= 0 && ibin < nbins) {
1104 l_int32 i, j, ns, nd, index, count, val;
1105 l_float32 start, oldsize;
1109 return (
NUMA *)ERROR_PTR(
"nas not defined", __func__, NULL);
1111 return (
NUMA *)ERROR_PTR(
"newsize must be > 1", __func__, NULL);
1113 return (
NUMA *)ERROR_PTR(
"no bins in nas", __func__, NULL);
1115 nd = (ns + newsize - 1) / newsize;
1117 return (
NUMA *)ERROR_PTR(
"nad not made", __func__, NULL);
1121 for (i = 0; i < nd; i++) {
1123 index = i * newsize;
1124 for (j = 0; j < newsize; j++) {
1151 l_float32 sum, factor, fval;
1155 return (
NUMA *)ERROR_PTR(
"nas not defined", __func__, NULL);
1157 return (
NUMA *)ERROR_PTR(
"tsum must be > 0.0", __func__, NULL);
1159 return (
NUMA *)ERROR_PTR(
"no bins in nas", __func__, NULL);
1162 factor = tsum / sum;
1165 return (
NUMA *)ERROR_PTR(
"nad not made", __func__, NULL);
1168 for (i = 0; i < ns; i++) {
1232 l_float32 *pvariance,
1239 l_float32 minval, maxval, fval, mean, sum;
1242 if (pmin) *pmin = 0.0;
1243 if (pmax) *pmax = 0.0;
1244 if (pmean) *pmean = 0.0;
1245 if (pvariance) *pvariance = 0.0;
1246 if (pmedian) *pmedian = 0.0;
1247 if (prval) *prval = 0.0;
1248 if (phisto) *phisto = NULL;
1250 return ERROR_INT(
"na not defined", __func__, 1);
1252 return ERROR_INT(
"numa is empty", __func__, 1);
1256 if (pmin) *pmin = minval;
1257 if (pmax) *pmax = maxval;
1258 if (pmean || pvariance) {
1260 for (i = 0; i < n; i++) {
1264 mean = sum / (l_float32)n;
1265 if (pmean) *pmean = mean;
1269 for (i = 0; i < n; i++) {
1273 *pvariance = sum / (l_float32)n - mean * mean;
1276 if (!pmedian && !prval && !phisto)
1320 l_float32 *pxmedian,
1322 l_float32 *pxvariance)
1324 if (pxmean) *pxmean = 0.0;
1325 if (pxmedian) *pxmedian = 0.0;
1326 if (pxmode) *pxmode = 0.0;
1327 if (pxvariance) *pxvariance = 0.0;
1329 return ERROR_INT(
"nahisto not defined", __func__, 1);
1332 pxmean, pxmedian, pxmode,
1369 l_float32 *pxmedian,
1371 l_float32 *pxvariance)
1374 l_float32 sum, sumval, halfsum, moment, var, x, y, ymax;
1376 if (pxmean) *pxmean = 0.0;
1377 if (pxmedian) *pxmedian = 0.0;
1378 if (pxmode) *pxmode = 0.0;
1379 if (pxvariance) *pxvariance = 0.0;
1381 return ERROR_INT(
"nahisto not defined", __func__, 1);
1382 if (!pxmean && !pxmedian && !pxmode && !pxvariance)
1383 return ERROR_INT(
"nothing to compute", __func__, 1);
1386 ifirst = L_MAX(0, ifirst);
1387 if (ilast < 0) ilast = n - 1;
1389 return ERROR_INT(
"invalid ifirst", __func__, 1);
1391 L_WARNING(
"ilast = %d is beyond max index = %d; adjusting\n",
1392 __func__, ilast, n - 1);
1396 return ERROR_INT(
"ifirst > ilast", __func__, 1);
1397 for (sum = 0.0, moment = 0.0, var = 0.0, i = ifirst; i <= ilast ; i++) {
1398 x = startx + i * deltax;
1405 L_INFO(
"sum is 0\n", __func__);
1410 *pxmean = moment / sum;
1412 *pxvariance = var / sum - moment * moment / (sum * sum);
1415 halfsum = sum / 2.0;
1416 for (sumval = 0.0, i = ifirst; i <= ilast; i++) {
1419 if (sumval >= halfsum) {
1420 *pxmedian = startx + i * deltax;
1429 for (i = ifirst; i <= ilast; i++) {
1436 *pxmode = startx + imax * deltax;
1463 l_float32 sum, fval;
1466 if (pnax) *pnax = NULL;
1468 return ERROR_INT(
"&nay not defined", __func__, 1);
1471 return ERROR_INT(
"nasy not defined", __func__, 1);
1473 return ERROR_INT(
"no bins in nas", __func__, 1);
1481 for (i = 0; i < n; i++) {
1490 startx, startx + n * deltax, npts,
1525 l_int32 i, ibinval, n;
1526 l_float32 startval, binsize, binval, maxval, fractval, total, sum, val;
1529 return ERROR_INT(
"prank not defined", __func__, 1);
1532 return ERROR_INT(
"na not defined", __func__, 1);
1535 if (rval < startval)
1537 maxval = startval + n * binsize;
1538 if (rval > maxval) {
1543 binval = (rval - startval) / binsize;
1544 ibinval = (l_int32)binval;
1549 fractval = binval - (l_float32)ibinval;
1552 for (i = 0; i < ibinval; i++) {
1557 sum += fractval * val;
1559 *prank = sum / total;
1595 l_float32 startval, binsize, rankcount, total, sum, fract, val;
1598 return ERROR_INT(
"prval not defined", __func__, 1);
1601 return ERROR_INT(
"na not defined", __func__, 1);
1603 L_WARNING(
"rank < 0; setting to 0.0\n", __func__);
1607 L_WARNING(
"rank > 1.0; setting to 1.0\n", __func__);
1614 rankcount = rank * total;
1616 for (i = 0; i < n; i++) {
1618 if (sum + val >= rankcount)
1625 fract = (rankcount - sum) / val;
1629 *prval = startval + binsize * ((l_float32)i + fract);
1666 l_int32 i, ntot, bincount, binindex, binsize;
1667 l_float32 sum, val, ave;
1670 return ERROR_INT(
"&nabinval not defined", __func__, 1);
1673 return ERROR_INT(
"na not defined", __func__, 1);
1675 return ERROR_INT(
"nbins must be > 1", __func__, 1);
1680 return ERROR_INT(
"naeach not made", __func__, 1);
1688 for (i = 0; i < ntot; i++) {
1692 if (bincount == binsize) {
1693 ave = sum / binsize;
1698 if (binindex == nbins)
break;
1702 *pnabinval = nabinval;
1742 l_int32 i, j, nxvals, occup, count, bincount, binindex, binsize;
1743 l_float32 sum, ave, ntot;
1745 if (pnarank) *pnarank = NULL;
1747 return ERROR_INT(
"&nabinval not defined", __func__, 1);
1750 return ERROR_INT(
"na not defined", __func__, 1);
1752 return ERROR_INT(
"nbins must be > 1", __func__, 1);
1756 occup = ntot / nxvals;
1757 if (occup < 1) L_INFO(
"average occupancy %d < 1\n", __func__, occup);
1761 return ERROR_INT(
"naeach not made", __func__, 1);
1769 for (i = 0; i < nxvals; i++) {
1771 for (j = 0; j < count; j++) {
1774 if (bincount == binsize) {
1775 ave = sum / binsize;
1780 if (binindex == nbins)
break;
1784 if (binindex == nbins)
break;
1786 *pnabinval = nabinval;
1787 if (binindex != nbins)
1788 L_ERROR(
"binindex = %d != nbins = %d\n", __func__, binindex, nbins);
1826 l_int32 maxbins, type;
1827 l_float32 maxval, delx;
1830 return ERROR_INT(
"&pnam not defined", __func__, 1);
1833 return ERROR_INT(
"na not defined", __func__, 1);
1835 return ERROR_INT(
"na is empty", __func__, 1);
1837 return ERROR_INT(
"nbins must be > 1", __func__, 1);
1847 L_INFO(
"sort the array: input size = %d\n", __func__,
numaGetCount(na));
1858 L_INFO(
"use a histogram: input size = %d\n", __func__,
numaGetCount(na));
1860 maxbins = L_MIN(100002, (l_int32)maxval + 2);
1867 L_WARNING(
"scale change: delx = %6.2f\n", __func__, delx);
1893 l_int32 i, start, end;
1897 return (
NUMA *)ERROR_PTR(
"ntotal <= 0", __func__, NULL);
1899 return (
NUMA *)ERROR_PTR(
"nbins <= 0", __func__, NULL);
1902 return (
NUMA *)ERROR_PTR(
"naeach not made", __func__, NULL);
1904 for (i = 0; i < nbins; i++) {
1905 end = ntotal * (i + 1) / nbins;
1967 l_float32 scorefract,
1968 l_int32 *psplitindex,
1975 l_int32 i, n, bestsplit, minrange, maxrange, maxindex;
1976 l_float32 ave1, ave2, ave1prev, ave2prev;
1977 l_float32 num1, num2, num1prev, num2prev;
1978 l_float32 val, minval, sum, fract1;
1979 l_float32 norm, score, minscore, maxscore;
1980 NUMA *nascore, *naave1, *naave2, *nanum1, *nanum2;
1982 if (psplitindex) *psplitindex = 0;
1983 if (pave1) *pave1 = 0.0;
1984 if (pave2) *pave2 = 0.0;
1985 if (pnum1) *pnum1 = 0.0;
1986 if (pnum2) *pnum2 = 0.0;
1987 if (pnascore) *pnascore = NULL;
1989 return ERROR_INT(
"na not defined", __func__, 1);
1993 return ERROR_INT(
"n = 1 in histogram", __func__, 1);
1996 return ERROR_INT(
"sum <= 0.0", __func__, 1);
1997 norm = 4.0 / ((l_float32)(n - 1) * (n - 1));
2008 return ERROR_INT(
"nascore not made", __func__, 1);
2014 for (i = 0; i < n; i++) {
2016 num1 = num1prev + val;
2020 ave1 = (num1prev * ave1prev + i * val) / num1;
2021 num2 = num2prev - val;
2025 ave2 = (num2prev * ave2prev - i * val) / num2;
2026 fract1 = num1 / sum;
2027 score = norm * (fract1 * (1 - fract1)) * (ave2 - ave1) * (ave2 - ave1);
2033 if (score > maxscore) {
2046 minscore = (1. - scorefract) * maxscore;
2047 for (i = maxindex - 1; i >= 0; i--) {
2053 for (i = maxindex + 1; i < n; i++) {
2060 bestsplit = minrange;
2061 for (i = minrange + 1; i <= maxrange; i++) {
2072 bestsplit = L_MIN(255, bestsplit + 1);
2074 if (psplitindex) *psplitindex = bestsplit;
2081 lept_stderr(
"minrange = %d, maxrange = %d\n", minrange, maxrange);
2084 "Score for split distribution");
2085 *pnascore = nascore;
2132 NUMA *na1, *na2, *nad;
2135 return ERROR_INT(
"&nad not defined", __func__, 1);
2138 return ERROR_INT(
"na1 and na2 not both defined", __func__, 1);
2141 return ERROR_INT(
"naa1 and naa2 numa counts differ", __func__, 1);
2144 return ERROR_INT(
"naa1 and naa2 number counts differ", __func__, 1);
2146 return ERROR_INT(
"na sizes must be 256", __func__, 1);
2150 for (i = 0; i < n; i++) {
2195 l_float32 sum1, sum2, diff, total;
2196 l_float32 *array1, *array3;
2200 return ERROR_INT(
"&dist not defined", __func__, 1);
2203 return ERROR_INT(
"na1 and na2 not both defined", __func__, 1);
2206 return ERROR_INT(
"na1 and na2 have different size", __func__, 1);
2211 norm = (L_ABS(sum1 - sum2) < 0.00001 * L_ABS(sum1)) ? 1 : 0;
2221 for (i = 1; i < n; i++) {
2222 diff = array1[i - 1] - array3[i - 1];
2224 total += L_ABS(diff);
2226 *pdist = total / sum1;
2286 l_int32 i, j, n, nn;
2288 l_float32 mean, var, rvar;
2289 NUMA *na1, *na2, *na3, *na4;
2291 if (pnam) *pnam = NULL;
2292 if (pnams) *pnams = NULL;
2293 if (pnav) *pnav = NULL;
2294 if (pnarv) *pnarv = NULL;
2295 if (!pnam && !pnams && !pnav && !pnarv)
2296 return ERROR_INT(
"nothing requested", __func__, 1);
2298 return ERROR_INT(
"naa not defined", __func__, 1);
2300 for (i = 0; i < n; i++) {
2303 L_ERROR(
"%d numbers in numa[%d]\n", __func__, nn, i);
2315 arrays = (l_float32 **)LEPT_CALLOC(n,
sizeof(l_float32 *));
2316 for (i = 0; i < n; i++) {
2327 for (j = 0; j < 256; j++) {
2329 for (i = 0; i < n; i++) {
2340 for (i = 0; i < n; i++)
2341 LEPT_FREE(arrays[i]);
2372 l_int32 i, k, n, maxloc, lloc, rloc;
2373 l_float32 fmaxval, sum, total, newtotal, val, lastval;
2374 l_float32 peakfract;
2378 return (
NUMA *)ERROR_PTR(
"nas not defined", __func__, NULL);
2384 return (
NUMA *)ERROR_PTR(
"na not made", __func__, NULL);
2385 if ((napeak =
numaCreate(4 * nmax)) == NULL) {
2387 return (
NUMA *)ERROR_PTR(
"napeak not made", __func__, NULL);
2390 for (k = 0; k < nmax; k++) {
2392 if (newtotal == 0.0)
2398 for (i = maxloc - 1; i >= 0; --i) {
2404 if (val > fract1 * fmaxval) {
2409 if (lastval - val > fract2 * lastval) {
2419 for (i = maxloc + 1; i < n; ++i) {
2425 if (val > fract1 * fmaxval) {
2430 if (lastval - val > fract2 * lastval) {
2438 peakfract = sum / total;
2444 for (i = lloc; i <= rloc; i++)
2487 l_int32 i, n, found, loc, direction;
2488 l_float32 startval, val, maxval, minval;
2491 if (pnav) *pnav = NULL;
2493 return (
NUMA *)ERROR_PTR(
"nas not defined", __func__, NULL);
2495 return (
NUMA *)ERROR_PTR(
"delta < 0", __func__, NULL);
2510 for (i = 1; i < n; i++) {
2512 if (L_ABS(val - startval) >= delta) {
2522 if (val > startval) {
2533 for (i = i + 1; i < n; i++) {
2535 if (direction == 1 && val > maxval ) {
2538 }
else if (direction == -1 && val < minval ) {
2541 }
else if (direction == 1 && (maxval - val >= delta)) {
2547 }
else if (direction == -1 && (val - minval >= delta)) {
2594 l_int32 i, n, start, index, minloc, found;
2595 l_float32 val, pval, jval, minval, maxval, sum, partsum;
2598 if (pfract) *pfract = 0.0;
2600 return ERROR_INT(
"&thresh not defined", __func__, 1);
2603 return ERROR_INT(
"na not defined", __func__, 1);
2604 if (skip <= 0) skip = 20;
2609 if (minval == maxval)
2610 return ERROR_INT(
"all array values are the same", __func__, 1);
2615 L_WARNING(
"array size %d < 256\n", __func__, n);
2618 for (i = 1; i < n; i++) {
2620 index = L_MIN(i + skip, n - 1);
2622 if (val < pval && jval < pval)
2628 return ERROR_INT(
"top of first peak not found", __func__, 1);
2634 for (i = start + 1; i < n; i++) {
2639 index = L_MIN(i + skip, n - 1);
2651 return ERROR_INT(
"no minimum found", __func__, 1);
2656 for (i = index - 1; i > index - skip; i--) {
2657 if (fa[i] < minval) {
2664 if (minloc > n - 10)
2665 return ERROR_INT(
"minimum at end of array; invalid", __func__, 1);
2673 *pfract = partsum / sum;
2702 l_float32 minreversal,
2706 l_int32 i, n, nr, ival, binvals;
2708 l_float32 fval, delx, len;
2712 if (prd) *prd = 0.0;
2714 return ERROR_INT(
"neither &nr nor &rd are defined", __func__, 1);
2716 return ERROR_INT(
"nas not defined", __func__, 1);
2718 L_INFO(
"nas is empty\n", __func__);
2721 if (minreversal < 0.0)
2722 return ERROR_INT(
"minreversal < 0", __func__, 1);
2726 for (i = 0; i < n; i++) {
2728 if (fval != 0.0 && fval != 1.0) {
2736 if (minreversal > 1.0) {
2737 L_WARNING(
"binary values but minreversal > 1\n", __func__);
2741 for (i = 1; i < n; i++) {
2742 if (ia[i] != ival) {
2758 *prd = (l_float32)nr / len;
2797 l_float32 estthresh,
2798 l_float32 *pbestthresh)
2800 l_int32 i, inrun, istart, iend, maxstart, maxend, runlen, maxrunlen;
2801 l_int32 val, maxval, nmax, count;
2802 l_float32 thresh, fmaxval, fmodeval;
2806 return ERROR_INT(
"&bestthresh not defined", __func__, 1);
2809 return ERROR_INT(
"nay not defined", __func__, 1);
2811 L_WARNING(
"nay count < 2; no threshold crossing\n", __func__);
2817 for (i = 0; i < 41; i++) {
2818 thresh = estthresh - 80.0 + 4.0 * i;
2827 maxval = (l_int32)fmaxval;
2829 for (i = 0; i < 41; i++) {
2836 if (count > nmax && fmodeval > 0.5 * fmaxval)
2837 maxval = (l_int32)fmodeval;
2842 maxrunlen = 0, maxstart = 0, maxend = 0;
2843 for (i = 0; i < 41; i++) {
2845 if (val == maxval) {
2852 if (inrun && (val != maxval)) {
2854 runlen = iend - istart + 1;
2856 if (runlen > maxrunlen) {
2864 runlen = i - istart;
2865 if (runlen > maxrunlen) {
2872 *pbestthresh = estthresh - 80.0 + 2.0 * (l_float32)(maxstart + maxend);
2875 lept_stderr(
"\nCrossings attain a maximum at %d thresholds, between:\n"
2876 " thresh[%d] = %5.1f and thresh[%d] = %5.1f\n",
2877 nmax, maxstart, estthresh - 80.0 + 4.0 * maxstart,
2878 maxend, estthresh - 80.0 + 4.0 * maxend);
2879 lept_stderr(
"The best choice: %5.1f\n", *pbestthresh);
2880 lept_stderr(
"Number of crossings at the 41 thresholds:");
2909 l_float32 startx, delx;
2910 l_float32 xval1, xval2, yval1, yval2, delta1, delta2, crossval, fract;
2914 return (
NUMA *)ERROR_PTR(
"nay not defined", __func__, NULL);
2918 return (
NUMA *)ERROR_PTR(
"nax and nay sizes differ", __func__, NULL);
2921 if (n < 2)
return nad;
2928 for (i = 1; i < n; i++) {
2933 xval2 = startx + i * delx;
2934 delta1 = yval1 - thresh;
2935 delta2 = yval2 - thresh;
2936 if (delta1 == 0.0) {
2938 }
else if (delta2 == 0.0) {
2940 }
else if (delta1 * delta2 < 0.0) {
2941 fract = L_ABS(delta1) / L_ABS(yval1 - yval2);
2942 crossval = xval1 + fract * (xval2 - xval1);
2972 l_int32 i, j, n, np, previndex, curindex;
2973 l_float32 startx, delx;
2974 l_float32 xval1, xval2, yval1, yval2, delta1, delta2;
2975 l_float32 prevval, curval, thresh, crossval, fract;
2979 return (
NUMA *)ERROR_PTR(
"nay not defined", __func__, NULL);
2983 return (
NUMA *)ERROR_PTR(
"nax and nay sizes differ", __func__, NULL);
2991 L_INFO(
"Number of crossings: %d\n", __func__, np);
2998 for (i = 0; i < np; i++) {
3001 thresh = (prevval + curval) / 2.0;
3005 xval1 = startx + previndex * delx;
3007 for (j = previndex + 1; j <= curindex; j++) {
3011 xval2 = startx + j * delx;
3013 delta1 = yval1 - thresh;
3014 delta2 = yval2 - thresh;
3015 if (delta1 == 0.0) {
3018 }
else if (delta2 == 0.0) {
3021 }
else if (delta1 * delta2 < 0.0) {
3022 fract = L_ABS(delta1) / L_ABS(yval1 - yval2);
3023 crossval = xval1 + fract * (xval2 - xval1);
3030 previndex = curindex;
3078 l_float32 relweight,
3083 l_float32 *pbestwidth,
3084 l_float32 *pbestshift,
3085 l_float32 *pbestscore)
3088 l_float32 delwidth, delshift, width, shift, score;
3089 l_float32 bestwidth, bestshift, bestscore;
3091 if (pbestscore) *pbestscore = 0.0;
3092 if (pbestwidth) *pbestwidth = 0.0;
3093 if (pbestshift) *pbestshift = 0.0;
3094 if (!pbestwidth || !pbestshift)
3095 return ERROR_INT(
"&bestwidth and &bestshift not defined", __func__, 1);
3097 return ERROR_INT(
"nas not defined", __func__, 1);
3099 bestscore = bestwidth = bestshift = 0.0;
3100 delwidth = (maxwidth - minwidth) / (nwidth - 1.0);
3101 for (i = 0; i < nwidth; i++) {
3102 width = minwidth + delwidth * i;
3103 delshift = width / (l_float32)(nshift);
3104 for (j = 0; j < nshift; j++) {
3105 shift = j * delshift;
3107 if (score > bestscore) {
3112 lept_stderr(
"width = %7.3f, shift = %7.3f, score = %7.3f\n",
3113 width, shift, score);
3119 *pbestwidth = bestwidth;
3120 *pbestshift = bestshift;
3122 *pbestscore = bestscore;
3163 l_float32 relweight,
3166 l_int32 i, n, nsamp, index;
3167 l_float32 score, weight, val;
3170 return ERROR_INT(
"&score not defined", __func__, 1);
3173 return ERROR_INT(
"nas not defined", __func__, 1);
3175 return ERROR_INT(
"nas size too small", __func__, 1);
3178 nsamp = (l_int32)((n - shift) / width);
3179 for (i = 0; i < nsamp; i++) {
3180 index = (l_int32)(shift + i * width);
3181 weight = (i % 2) ? 1.0 : -1.0 * relweight;
3183 score += weight * val;
3186 *pscore = 2.0 * width * score / (l_float32)n;
3220 l_int32 i, nsets, val;
3224 first = L_MAX(0, first);
3226 return (
NUMA *)ERROR_PTR(
"last < first!", __func__, NULL);
3228 return (
NUMA *)ERROR_PTR(
"nmax < 1!", __func__, NULL);
3230 nsets = L_MIN(nmax, last - first + 1);
3234 return (
NUMA *)ERROR_PTR(
"nsets == 0", __func__, NULL);
3241 delta = (l_float32)(last - first) / (nsets - 1);
3243 delta = (l_float32)(last - first - 1) / (nsets - 1);
3247 for (i = 0; i < nsets; i++) {
3248 val = (l_int32)(first + i * delta + 0.5);
l_ok gplotSimple1(NUMA *na, l_int32 outformat, const char *outroot, const char *title)
gplotSimple1()
l_ok numaAddNumber(NUMA *na, l_float32 val)
numaAddNumber()
l_int32 numaaGetNumberCount(NUMAA *naa)
numaaGetNumberCount()
l_ok numaGetFValue(NUMA *na, l_int32 index, l_float32 *pval)
numaGetFValue()
l_ok numaWriteStderr(NUMA *na)
numaWriteStderr()
NUMA * numaCreate(l_int32 n)
numaCreate()
l_int32 numaaGetCount(NUMAA *naa)
numaaGetCount()
NUMA * numaaGetNuma(NUMAA *naa, l_int32 index, l_int32 accessflag)
numaaGetNuma()
l_int32 numaaGetNumaCount(NUMAA *naa, l_int32 index)
numaaGetNumaCount()
l_ok numaSetCount(NUMA *na, l_int32 newcount)
numaSetCount()
void numaDestroy(NUMA **pna)
numaDestroy()
l_ok numaSetValue(NUMA *na, l_int32 index, l_float32 val)
numaSetValue()
l_int32 numaGetCount(NUMA *na)
numaGetCount()
l_ok numaGetIValue(NUMA *na, l_int32 index, l_int32 *pival)
numaGetIValue()
NUMA * numaCopy(NUMA *na)
numaCopy()
l_int32 * numaGetIArray(NUMA *na)
numaGetIArray()
l_ok numaGetParameters(NUMA *na, l_float32 *pstartx, l_float32 *pdelx)
numaGetParameters()
l_ok numaCopyParameters(NUMA *nad, NUMA *nas)
numaCopyParameters()
l_ok numaSetParameters(NUMA *na, l_float32 startx, l_float32 delx)
numaSetParameters()
l_float32 * numaGetFArray(NUMA *na, l_int32 copyflag)
numaGetFArray()
NUMA * numaSort(NUMA *naout, NUMA *nain, l_int32 sortorder)
numaSort()
l_ok numaGetMode(NUMA *na, l_float32 *pval, l_int32 *pcount)
numaGetMode()
NUMA * numaAddBorder(NUMA *nas, l_int32 left, l_int32 right, l_float32 val)
numaAddBorder()
l_ok numaHasOnlyIntegers(NUMA *na, l_int32 *pallints)
numaHasOnlyIntegers()
l_ok numaGetMedian(NUMA *na, l_float32 *pval)
numaGetMedian()
l_ok numaGetSumOnInterval(NUMA *na, l_int32 first, l_int32 last, l_float32 *psum)
numaGetSumOnInterval()
NUMA * numaRemoveBorder(NUMA *nas, l_int32 left, l_int32 right)
numaRemoveBorder()
l_ok numaInterpolateEqxInterval(l_float32 startx, l_float32 deltax, NUMA *nasy, l_int32 type, l_float32 x0, l_float32 x1, l_int32 npts, NUMA **pnax, NUMA **pnay)
numaInterpolateEqxInterval()
l_ok numaGetMin(NUMA *na, l_float32 *pminval, l_int32 *piminloc)
numaGetMin()
NUMA * numaMakeConstant(l_float32 val, l_int32 size)
numaMakeConstant()
NUMA * numaClipToInterval(NUMA *nas, l_int32 first, l_int32 last)
numaClipToInterval()
l_ok numaGetMax(NUMA *na, l_float32 *pmaxval, l_int32 *pimaxloc)
numaGetMax()
NUMA * numaGetPartialSums(NUMA *na)
numaGetPartialSums()
l_int32 numaChooseSortType(NUMA *nas)
numaChooseSortType()
NUMA * numaAddSpecifiedBorder(NUMA *nas, l_int32 left, l_int32 right, l_int32 type)
numaAddSpecifiedBorder()
l_ok numaGetSum(NUMA *na, l_float32 *psum)
numaGetSum()
NUMA * numaWindowedMeanSquare(NUMA *nas, l_int32 wc)
numaWindowedMeanSquare()
l_ok numaGetRankBinValues(NUMA *na, l_int32 nbins, NUMA **pnam)
numaGetRankBinValues()
l_ok numaCountReversals(NUMA *nas, l_float32 minreversal, l_int32 *pnr, l_float32 *prd)
numaCountReversals()
NUMA * numaMakeHistogram(NUMA *na, l_int32 maxbins, l_int32 *pbinsize, l_int32 *pbinstart)
numaMakeHistogram()
NUMA * numaFindExtrema(NUMA *nas, l_float32 delta, NUMA **pnav)
numaFindExtrema()
NUMA * numaDilate(NUMA *nas, l_int32 size)
numaDilate()
NUMA * genConstrainedNumaInRange(l_int32 first, l_int32 last, l_int32 nmax, l_int32 use_pairs)
genConstrainedNumaInRange()
l_ok numaEvalHaarSum(NUMA *nas, l_float32 width, l_float32 shift, l_float32 relweight, l_float32 *pscore)
numaEvalHaarSum()
NUMA * numaConvertToInt(NUMA *nas)
numaConvertToInt()
l_ok numaEarthMoverDistance(NUMA *na1, NUMA *na2, l_float32 *pdist)
numaEarthMoverDistance()
l_ok numaSelectCrossingThreshold(NUMA *nax, NUMA *nay, l_float32 estthresh, l_float32 *pbestthresh)
numaSelectCrossingThreshold()
l_ok numaWindowedStats(NUMA *nas, l_int32 wc, NUMA **pnam, NUMA **pnams, NUMA **pnav, NUMA **pnarv)
numaWindowedStats()
NUMA * numaFindPeaks(NUMA *nas, l_int32 nmax, l_float32 fract1, l_float32 fract2)
numaFindPeaks()
NUMA * numaWindowedMean(NUMA *nas, l_int32 wc)
numaWindowedMean()
l_ok numaWindowedVariance(NUMA *nam, NUMA *nams, NUMA **pnav, NUMA **pnarv)
numaWindowedVariance()
l_ok numaEvalBestHaarParameters(NUMA *nas, l_float32 relweight, l_int32 nwidth, l_int32 nshift, l_float32 minwidth, l_float32 maxwidth, l_float32 *pbestwidth, l_float32 *pbestshift, l_float32 *pbestscore)
numaEvalBestHaarParameters()
l_ok numaDiscretizeHistoInBins(NUMA *na, l_int32 nbins, NUMA **pnabinval, NUMA **pnarank)
numaDiscretizeHistoInBins()
l_ok numaGetHistogramStatsOnInterval(NUMA *nahisto, l_float32 startx, l_float32 deltax, l_int32 ifirst, l_int32 ilast, l_float32 *pxmean, l_float32 *pxmedian, l_float32 *pxmode, l_float32 *pxvariance)
numaGetHistogramStatsOnInterval()
NUMA * numaMakeHistogramClipped(NUMA *na, l_float32 binsize, l_float32 maxsize)
numaMakeHistogramClipped()
NUMA * numaNormalizeHistogram(NUMA *nas, l_float32 tsum)
numaNormalizeHistogram()
NUMA * numaErode(NUMA *nas, l_int32 size)
numaErode()
l_ok numaGetHistogramStats(NUMA *nahisto, l_float32 startx, l_float32 deltax, l_float32 *pxmean, l_float32 *pxmedian, l_float32 *pxmode, l_float32 *pxvariance)
numaGetHistogramStats()
l_ok grayHistogramsToEMD(NUMAA *naa1, NUMAA *naa2, NUMA **pnad)
grayHistogramsToEMD()
NUMA * numaTransform(NUMA *nas, l_float32 shift, l_float32 scale)
numaTransform()
NUMA * numaWindowedMedian(NUMA *nas, l_int32 halfwin)
numaWindowedMedian()
l_ok numaMakeRankFromHistogram(l_float32 startx, l_float32 deltax, NUMA *nasy, l_int32 npts, NUMA **pnax, NUMA **pnay)
numaMakeRankFromHistogram()
l_ok numaSplitDistribution(NUMA *na, l_float32 scorefract, l_int32 *psplitindex, l_float32 *pave1, l_float32 *pave2, l_float32 *pnum1, l_float32 *pnum2, NUMA **pnascore)
numaSplitDistribution()
l_ok numaHistogramGetRankFromVal(NUMA *na, l_float32 rval, l_float32 *prank)
numaHistogramGetRankFromVal()
l_ok numaHistogramGetValFromRank(NUMA *na, l_float32 rank, l_float32 *prval)
numaHistogramGetValFromRank()
l_ok numaFindLocForThreshold(NUMA *na, l_int32 skip, l_int32 *pthresh, l_float32 *pfract)
numaFindLocForThreshold()
l_ok grayInterHistogramStats(NUMAA *naa, l_int32 wc, NUMA **pnam, NUMA **pnams, NUMA **pnav, NUMA **pnarv)
grayInterHistogramStats()
l_ok numaSimpleStats(NUMA *na, l_int32 first, l_int32 last, l_float32 *pmean, l_float32 *pvar, l_float32 *prvar)
numaSimpleStats()
NUMA * numaGetUniformBinSizes(l_int32 ntotal, l_int32 nbins)
numaGetUniformBinSizes()
NUMA * numaCrossingsByThreshold(NUMA *nax, NUMA *nay, l_float32 thresh)
numaCrossingsByThreshold()
NUMA * numaRebinHistogram(NUMA *nas, l_int32 newsize)
numaRebinHistogram()
NUMA * numaMakeHistogramAuto(NUMA *na, l_int32 maxbins)
numaMakeHistogramAuto()
NUMA * numaClose(NUMA *nas, l_int32 size)
numaClose()
l_ok numaDiscretizeSortedInBins(NUMA *na, l_int32 nbins, NUMA **pnabinval)
numaDiscretizeSortedInBins()
NUMA * numaOpen(NUMA *nas, l_int32 size)
numaOpen()
l_ok numaGetStatsUsingHistogram(NUMA *na, l_int32 maxbins, l_float32 *pmin, l_float32 *pmax, l_float32 *pmean, l_float32 *pvariance, l_float32 *pmedian, l_float32 rank, l_float32 *prval, NUMA **phisto)
numaGetStatsUsingHistogram()
NUMA * numaCrossingsByPeaks(NUMA *nax, NUMA *nay, l_float32 delta)
numaCrossingsByPeaks()
void lept_stderr(const char *fmt,...)
lept_stderr()