Leptonica  1.82.0
Image processing and image analysis suite
numafunc1.c
Go to the documentation of this file.
1 /*====================================================================*
2  - Copyright (C) 2001 Leptonica. All rights reserved.
3  -
4  - Redistribution and use in source and binary forms, with or without
5  - modification, are permitted provided that the following conditions
6  - are met:
7  - 1. Redistributions of source code must retain the above copyright
8  - notice, this list of conditions and the following disclaimer.
9  - 2. Redistributions in binary form must reproduce the above
10  - copyright notice, this list of conditions and the following
11  - disclaimer in the documentation and/or other materials
12  - provided with the distribution.
13  -
14  - THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
15  - ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
16  - LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
17  - A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL ANY
18  - CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
19  - EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
20  - PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
21  - PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY
22  - OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
23  - NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
24  - SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
25  *====================================================================*/
26 
144 #ifdef HAVE_CONFIG_H
145 #include <config_auto.h>
146 #endif /* HAVE_CONFIG_H */
147 
148 #include <math.h>
149 #include "allheaders.h"
150 
151 /*----------------------------------------------------------------------*
152  * Arithmetic and logical ops on Numas *
153  *----------------------------------------------------------------------*/
172 NUMA *
174  NUMA *na1,
175  NUMA *na2,
176  l_int32 op)
177 {
178 l_int32 i, n;
179 l_float32 val1, val2;
180 
181  PROCNAME("numaArithOp");
182 
183  if (!na1 || !na2)
184  return (NUMA *)ERROR_PTR("na1, na2 not both defined", procName, nad);
185  n = numaGetCount(na1);
186  if (n != numaGetCount(na2))
187  return (NUMA *)ERROR_PTR("na1, na2 sizes differ", procName, nad);
188  if (nad && nad != na1)
189  return (NUMA *)ERROR_PTR("nad defined but not in-place", procName, nad);
190  if (op != L_ARITH_ADD && op != L_ARITH_SUBTRACT &&
191  op != L_ARITH_MULTIPLY && op != L_ARITH_DIVIDE)
192  return (NUMA *)ERROR_PTR("invalid op", procName, nad);
193  if (op == L_ARITH_DIVIDE) {
194  for (i = 0; i < n; i++) {
195  numaGetFValue(na2, i, &val2);
196  if (val2 == 0.0)
197  return (NUMA *)ERROR_PTR("na2 has 0 element", procName, nad);
198  }
199  }
200 
201  /* If nad is not identical to na1, make it an identical copy */
202  if (!nad)
203  nad = numaCopy(na1);
204 
205  for (i = 0; i < n; i++) {
206  numaGetFValue(nad, i, &val1);
207  numaGetFValue(na2, i, &val2);
208  switch (op) {
209  case L_ARITH_ADD:
210  numaSetValue(nad, i, val1 + val2);
211  break;
212  case L_ARITH_SUBTRACT:
213  numaSetValue(nad, i, val1 - val2);
214  break;
215  case L_ARITH_MULTIPLY:
216  numaSetValue(nad, i, val1 * val2);
217  break;
218  case L_ARITH_DIVIDE:
219  numaSetValue(nad, i, val1 / val2);
220  break;
221  default:
222  lept_stderr(" Unknown arith op: %d\n", op);
223  return nad;
224  }
225  }
226 
227  return nad;
228 }
229 
230 
252 NUMA *
254  NUMA *na1,
255  NUMA *na2,
256  l_int32 op)
257 {
258 l_int32 i, n, val1, val2, val;
259 
260  PROCNAME("numaLogicalOp");
261 
262  if (!na1 || !na2)
263  return (NUMA *)ERROR_PTR("na1, na2 not both defined", procName, nad);
264  n = numaGetCount(na1);
265  if (n != numaGetCount(na2))
266  return (NUMA *)ERROR_PTR("na1, na2 sizes differ", procName, nad);
267  if (nad && nad != na1)
268  return (NUMA *)ERROR_PTR("nad defined; not in-place", procName, nad);
269  if (op != L_UNION && op != L_INTERSECTION &&
270  op != L_SUBTRACTION && op != L_EXCLUSIVE_OR)
271  return (NUMA *)ERROR_PTR("invalid op", procName, nad);
272 
273  /* If nad is not identical to na1, make it an identical copy */
274  if (!nad)
275  nad = numaCopy(na1);
276 
277  for (i = 0; i < n; i++) {
278  numaGetIValue(nad, i, &val1);
279  numaGetIValue(na2, i, &val2);
280  val1 = (val1 == 0) ? 0 : 1;
281  val2 = (val2 == 0) ? 0 : 1;
282  switch (op) {
283  case L_UNION:
284  val = (val1 || val2) ? 1 : 0;
285  numaSetValue(nad, i, val);
286  break;
287  case L_INTERSECTION:
288  val = (val1 && val2) ? 1 : 0;
289  numaSetValue(nad, i, val);
290  break;
291  case L_SUBTRACTION:
292  val = (val1 && !val2) ? 1 : 0;
293  numaSetValue(nad, i, val);
294  break;
295  case L_EXCLUSIVE_OR:
296  val = (val1 != val2) ? 1 : 0;
297  numaSetValue(nad, i, val);
298  break;
299  default:
300  lept_stderr(" Unknown logical op: %d\n", op);
301  return nad;
302  }
303  }
304 
305  return nad;
306 }
307 
308 
325 NUMA *
327  NUMA *nas)
328 {
329 l_int32 i, n, val;
330 
331  PROCNAME("numaInvert");
332 
333  if (!nas)
334  return (NUMA *)ERROR_PTR("nas not defined", procName, nad);
335  if (nad && nad != nas)
336  return (NUMA *)ERROR_PTR("nad defined; not in-place", procName, nad);
337 
338  if (!nad)
339  nad = numaCopy(nas);
340  n = numaGetCount(nad);
341  for (i = 0; i < n; i++) {
342  numaGetIValue(nad, i, &val);
343  if (!val)
344  val = 1;
345  else
346  val = 0;
347  numaSetValue(nad, i, val);
348  }
349  return nad;
350 }
351 
352 
369 l_int32
371  NUMA *na2,
372  l_float32 maxdiff,
373  l_int32 *psimilar)
374 {
375 l_int32 i, n;
376 l_float32 val1, val2;
377 
378  PROCNAME("numaSimilar");
379 
380  if (!psimilar)
381  return ERROR_INT("&similar not defined", procName, 1);
382  *psimilar = 0;
383  if (!na1 || !na2)
384  return ERROR_INT("na1 and na2 not both defined", procName, 1);
385  maxdiff = L_ABS(maxdiff);
386 
387  n = numaGetCount(na1);
388  if (n != numaGetCount(na2)) return 0;
389 
390  for (i = 0; i < n; i++) {
391  numaGetFValue(na1, i, &val1);
392  numaGetFValue(na2, i, &val2);
393  if (L_ABS(val1 - val2) > maxdiff) return 0;
394  }
395 
396  *psimilar = 1;
397  return 0;
398 }
399 
400 
418 l_ok
420  l_int32 index,
421  l_float32 val)
422 {
423 l_int32 n;
424 
425  PROCNAME("numaAddToNumber");
426 
427  if (!na)
428  return ERROR_INT("na not defined", procName, 1);
429  if ((n = numaGetCount(na)) == 0)
430  return ERROR_INT("na is empty", procName, 1);
431  if (index < 0 || index >= n) {
432  L_ERROR("index %d not in [0,...,%d]\n", procName, index, n - 1);
433  return 1;
434  }
435 
436  na->array[index] += val;
437  return 0;
438 }
439 
440 
441 /*----------------------------------------------------------------------*
442  * Simple extractions *
443  *----------------------------------------------------------------------*/
452 l_ok
454  l_float32 *pminval,
455  l_int32 *piminloc)
456 {
457 l_int32 i, n, iminloc;
458 l_float32 val, minval;
459 
460  PROCNAME("numaGetMin");
461 
462  if (!pminval && !piminloc)
463  return ERROR_INT("nothing to do", procName, 1);
464  if (pminval) *pminval = 0.0;
465  if (piminloc) *piminloc = 0;
466  if (!na)
467  return ERROR_INT("na not defined", procName, 1);
468  if ((n = numaGetCount(na)) == 0)
469  return ERROR_INT("na is empty", procName, 1);
470 
471  minval = +1000000000.;
472  iminloc = 0;
473  for (i = 0; i < n; i++) {
474  numaGetFValue(na, i, &val);
475  if (val < minval) {
476  minval = val;
477  iminloc = i;
478  }
479  }
480 
481  if (pminval) *pminval = minval;
482  if (piminloc) *piminloc = iminloc;
483  return 0;
484 }
485 
486 
495 l_ok
497  l_float32 *pmaxval,
498  l_int32 *pimaxloc)
499 {
500 l_int32 i, n, imaxloc;
501 l_float32 val, maxval;
502 
503  PROCNAME("numaGetMax");
504 
505  if (!pmaxval && !pimaxloc)
506  return ERROR_INT("nothing to do", procName, 1);
507  if (pmaxval) *pmaxval = 0.0;
508  if (pimaxloc) *pimaxloc = 0;
509  if (!na)
510  return ERROR_INT("na not defined", procName, 1);
511  if ((n = numaGetCount(na)) == 0)
512  return ERROR_INT("na is empty", procName, 1);
513 
514  maxval = -1000000000.;
515  imaxloc = 0;
516  for (i = 0; i < n; i++) {
517  numaGetFValue(na, i, &val);
518  if (val > maxval) {
519  maxval = val;
520  imaxloc = i;
521  }
522  }
523 
524  if (pmaxval) *pmaxval = maxval;
525  if (pimaxloc) *pimaxloc = imaxloc;
526  return 0;
527 }
528 
529 
537 l_ok
539  l_float32 *psum)
540 {
541 l_int32 i, n;
542 l_float32 val, sum;
543 
544  PROCNAME("numaGetSum");
545 
546  if (!psum)
547  return ERROR_INT("&sum not defined", procName, 1);
548  *psum = 0;
549  if (!na)
550  return ERROR_INT("na not defined", procName, 1);
551 
552  if ((n = numaGetCount(na)) == 0)
553  return ERROR_INT("na is empty", procName, 1);
554  sum = 0.0;
555  for (i = 0; i < n; i++) {
556  numaGetFValue(na, i, &val);
557  sum += val;
558  }
559  *psum = sum;
560  return 0;
561 }
562 
563 
578 NUMA *
580 {
581 l_int32 i, n;
582 l_float32 val, sum;
583 NUMA *nasum;
584 
585  PROCNAME("numaGetPartialSums");
586 
587  if (!na)
588  return (NUMA *)ERROR_PTR("na not defined", procName, NULL);
589 
590  if ((n = numaGetCount(na)) == 0)
591  L_WARNING("na is empty\n", procName);
592  nasum = numaCreate(n);
593  sum = 0.0;
594  for (i = 0; i < n; i++) {
595  numaGetFValue(na, i, &val);
596  sum += val;
597  numaAddNumber(nasum, sum);
598  }
599  return nasum;
600 }
601 
602 
612 l_ok
614  l_int32 first,
615  l_int32 last,
616  l_float32 *psum)
617 {
618 l_int32 i, n;
619 l_float32 val, sum;
620 
621  PROCNAME("numaGetSumOnInterval");
622 
623  if (!psum)
624  return ERROR_INT("&sum not defined", procName, 1);
625  *psum = 0.0;
626  if (!na)
627  return ERROR_INT("na not defined", procName, 1);
628 
629  sum = 0.0;
630  if ((n = numaGetCount(na)) == 0)
631  return ERROR_INT("na is empty", procName, 1);
632  if (first < 0) first = 0;
633  if (first >= n || last < -1) /* not an error */
634  return 0;
635  if (last == -1)
636  last = n - 1;
637  last = L_MIN(last, n - 1);
638 
639  for (i = first; i <= last; i++) {
640  numaGetFValue(na, i, &val);
641  sum += val;
642  }
643  *psum = sum;
644  return 0;
645 }
646 
647 
655 l_ok
657  l_int32 *pallints)
658 {
659 l_int32 i, n;
660 l_float32 val;
661 
662  PROCNAME("numaHasOnlyIntegers");
663 
664  if (!pallints)
665  return ERROR_INT("&allints not defined", procName, 1);
666  *pallints = TRUE;
667  if (!na)
668  return ERROR_INT("na not defined", procName, 1);
669 
670  if ((n = numaGetCount(na)) == 0)
671  return ERROR_INT("na is empty", procName, 1);
672  for (i = 0; i < n; i ++) {
673  numaGetFValue(na, i, &val);
674  if (val != (l_int32)val) {
675  *pallints = FALSE;
676  return 0;
677  }
678  }
679  return 0;
680 }
681 
682 
690 l_ok
692  l_float32 *pave)
693 {
694 l_int32 n;
695 l_float32 sum;
696 
697  PROCNAME("numaGetMean");
698 
699  if (!pave)
700  return ERROR_INT("&ave not defined", procName, 1);
701  *pave = 0;
702  if (!na)
703  return ERROR_INT("na not defined", procName, 1);
704  if ((n = numaGetCount(na)) == 0)
705  return ERROR_INT("na is empty", procName, 1);
706 
707  numaGetSum(na, &sum);
708  *pave = sum / n;
709  return 0;
710 }
711 
719 l_ok
721  l_float32 *paveabs)
722 {
723 l_int32 n;
724 NUMA *na1;
725 
726  PROCNAME("numaGetMeanAbsval");
727 
728  if (!paveabs)
729  return ERROR_INT("&aveabs not defined", procName, 1);
730  *paveabs = 0;
731  if (!na)
732  return ERROR_INT("na not defined", procName, 1);
733  if ((n = numaGetCount(na)) == 0)
734  return ERROR_INT("na is empty", procName, 1);
735 
736  na1 = numaMakeAbsval(NULL, na);
737  numaGetMean(na1, paveabs);
738  numaDestroy(&na1);
739  return 0;
740 }
741 
742 
750 NUMA *
752  l_int32 subfactor)
753 {
754 l_int32 i, n;
755 l_float32 val;
756 NUMA *nad;
757 
758  PROCNAME("numaSubsample");
759 
760  if (!nas)
761  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
762  if (subfactor < 1)
763  return (NUMA *)ERROR_PTR("subfactor < 1", procName, NULL);
764 
765  nad = numaCreate(0);
766  if ((n = numaGetCount(nas)) == 0)
767  L_WARNING("nas is empty\n", procName);
768  for (i = 0; i < n; i++) {
769  if (i % subfactor != 0) continue;
770  numaGetFValue(nas, i, &val);
771  numaAddNumber(nad, val);
772  }
773 
774  return nad;
775 }
776 
777 
785 NUMA *
787 {
788 l_int32 i, n;
789 l_float32 prev, cur;
790 NUMA *nad;
791 
792  PROCNAME("numaMakeDelta");
793 
794  if (!nas)
795  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
796  if ((n = numaGetCount(nas)) < 2) {
797  L_WARNING("n < 2; returning empty numa\n", procName);
798  return numaCreate(1);
799  }
800 
801  nad = numaCreate(n - 1);
802  numaGetFValue(nas, 0, &prev);
803  for (i = 1; i < n; i++) {
804  numaGetFValue(nas, i, &cur);
805  numaAddNumber(nad, cur - prev);
806  prev = cur;
807  }
808  return nad;
809 }
810 
811 
820 NUMA *
821 numaMakeSequence(l_float32 startval,
822  l_float32 increment,
823  l_int32 size)
824 {
825 l_int32 i;
826 l_float32 val;
827 NUMA *na;
828 
829  PROCNAME("numaMakeSequence");
830 
831  if ((na = numaCreate(size)) == NULL)
832  return (NUMA *)ERROR_PTR("na not made", procName, NULL);
833 
834  for (i = 0; i < size; i++) {
835  val = startval + i * increment;
836  numaAddNumber(na, val);
837  }
838  return na;
839 }
840 
841 
850 NUMA *
851 numaMakeConstant(l_float32 val,
852  l_int32 size)
853 {
854  return numaMakeSequence(val, 0.0, size);
855 }
856 
857 
866 NUMA *
868  NUMA *nas)
869 {
870 l_int32 i, n;
871 l_float32 val;
872 
873  PROCNAME("numaMakeAbsval");
874 
875  if (!nas)
876  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
877  if (nad && nad != nas)
878  return (NUMA *)ERROR_PTR("nad and not in-place", procName, NULL);
879 
880  if (!nad)
881  nad = numaCopy(nas);
882  n = numaGetCount(nad);
883  for (i = 0; i < n; i++) {
884  val = nad->array[i];
885  nad->array[i] = L_ABS(val);
886  }
887 
888  return nad;
889 }
890 
891 
901 NUMA *
903  l_int32 left,
904  l_int32 right,
905  l_float32 val)
906 {
907 l_int32 i, n, len;
908 l_float32 startx, delx;
909 l_float32 *fas, *fad;
910 NUMA *nad;
911 
912  PROCNAME("numaAddBorder");
913 
914  if (!nas)
915  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
916  if (left < 0) left = 0;
917  if (right < 0) right = 0;
918  if (left == 0 && right == 0)
919  return numaCopy(nas);
920 
921  n = numaGetCount(nas);
922  len = n + left + right;
923  nad = numaMakeConstant(val, len);
924  numaGetParameters(nas, &startx, &delx);
925  numaSetParameters(nad, startx - delx * left, delx);
926  fas = numaGetFArray(nas, L_NOCOPY);
927  fad = numaGetFArray(nad, L_NOCOPY);
928  for (i = 0; i < n; i++)
929  fad[left + i] = fas[i];
930 
931  return nad;
932 }
933 
934 
944 NUMA *
946  l_int32 left,
947  l_int32 right,
948  l_int32 type)
949 {
950 l_int32 i, n;
951 l_float32 *fa;
952 NUMA *nad;
953 
954  PROCNAME("numaAddSpecifiedBorder");
955 
956  if (!nas)
957  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
958  if (left < 0) left = 0;
959  if (right < 0) right = 0;
960  if (left == 0 && right == 0)
961  return numaCopy(nas);
962  if (type != L_CONTINUED_BORDER && type != L_MIRRORED_BORDER)
963  return (NUMA *)ERROR_PTR("invalid type", procName, NULL);
964  n = numaGetCount(nas);
965  if (type == L_MIRRORED_BORDER && (left > n || right > n))
966  return (NUMA *)ERROR_PTR("border too large", procName, NULL);
967 
968  nad = numaAddBorder(nas, left, right, 0);
969  n = numaGetCount(nad);
970  fa = numaGetFArray(nad, L_NOCOPY);
971  if (type == L_CONTINUED_BORDER) {
972  for (i = 0; i < left; i++)
973  fa[i] = fa[left];
974  for (i = n - right; i < n; i++)
975  fa[i] = fa[n - right - 1];
976  } else { /* type == L_MIRRORED_BORDER */
977  for (i = 0; i < left; i++)
978  fa[i] = fa[2 * left - 1 - i];
979  for (i = 0; i < right; i++)
980  fa[n - right + i] = fa[n - right - i - 1];
981  }
982 
983  return nad;
984 }
985 
986 
995 NUMA *
997  l_int32 left,
998  l_int32 right)
999 {
1000 l_int32 i, n, len;
1001 l_float32 startx, delx;
1002 l_float32 *fas, *fad;
1003 NUMA *nad;
1004 
1005  PROCNAME("numaRemoveBorder");
1006 
1007  if (!nas)
1008  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
1009  if (left < 0) left = 0;
1010  if (right < 0) right = 0;
1011  if (left == 0 && right == 0)
1012  return numaCopy(nas);
1013 
1014  n = numaGetCount(nas);
1015  if ((len = n - left - right) < 0)
1016  return (NUMA *)ERROR_PTR("len < 0 after removal", procName, NULL);
1017  nad = numaMakeConstant(0, len);
1018  numaGetParameters(nas, &startx, &delx);
1019  numaSetParameters(nad, startx + delx * left, delx);
1020  fas = numaGetFArray(nas, L_NOCOPY);
1021  fad = numaGetFArray(nad, L_NOCOPY);
1022  for (i = 0; i < len; i++)
1023  fad[i] = fas[left + i];
1024 
1025  return nad;
1026 }
1027 
1028 
1036 l_ok
1038  l_int32 *pcount)
1039 {
1040 l_int32 n, i, val, count, inrun;
1041 
1042  PROCNAME("numaCountNonzeroRuns");
1043 
1044  if (!pcount)
1045  return ERROR_INT("&count not defined", procName, 1);
1046  *pcount = 0;
1047  if (!na)
1048  return ERROR_INT("na not defined", procName, 1);
1049  if ((n = numaGetCount(na)) == 0)
1050  return ERROR_INT("na is empty", procName, 1);
1051 
1052  count = 0;
1053  inrun = FALSE;
1054  for (i = 0; i < n; i++) {
1055  numaGetIValue(na, i, &val);
1056  if (!inrun && val > 0) {
1057  count++;
1058  inrun = TRUE;
1059  } else if (inrun && val == 0) {
1060  inrun = FALSE;
1061  }
1062  }
1063  *pcount = count;
1064  return 0;
1065 }
1066 
1067 
1077 l_ok
1079  l_float32 eps,
1080  l_int32 *pfirst,
1081  l_int32 *plast)
1082 {
1083 l_int32 n, i, found;
1084 l_float32 val;
1085 
1086  PROCNAME("numaGetNonzeroRange");
1087 
1088  if (pfirst) *pfirst = 0;
1089  if (plast) *plast = 0;
1090  if (!pfirst || !plast)
1091  return ERROR_INT("pfirst and plast not both defined", procName, 1);
1092  if (!na)
1093  return ERROR_INT("na not defined", procName, 1);
1094  if ((n = numaGetCount(na)) == 0)
1095  return ERROR_INT("na is empty", procName, 1);
1096 
1097  found = FALSE;
1098  for (i = 0; i < n; i++) {
1099  numaGetFValue(na, i, &val);
1100  if (val > eps) {
1101  found = TRUE;
1102  break;
1103  }
1104  }
1105  if (!found) {
1106  *pfirst = n - 1;
1107  *plast = 0;
1108  return 1;
1109  }
1110 
1111  *pfirst = i;
1112  for (i = n - 1; i >= 0; i--) {
1113  numaGetFValue(na, i, &val);
1114  if (val > eps)
1115  break;
1116  }
1117  *plast = i;
1118  return 0;
1119 }
1120 
1121 
1130 l_ok
1132  l_int32 type,
1133  l_int32 *pcount)
1134 {
1135 l_int32 n, i, count;
1136 l_float32 val;
1137 
1138  PROCNAME("numaGetCountRelativeToZero");
1139 
1140  if (!pcount)
1141  return ERROR_INT("&count not defined", procName, 1);
1142  *pcount = 0;
1143  if (!na)
1144  return ERROR_INT("na not defined", procName, 1);
1145  if ((n = numaGetCount(na)) == 0)
1146  return ERROR_INT("na is empty", procName, 1);
1147 
1148  for (i = 0, count = 0; i < n; i++) {
1149  numaGetFValue(na, i, &val);
1150  if (type == L_LESS_THAN_ZERO && val < 0.0)
1151  count++;
1152  else if (type == L_EQUAL_TO_ZERO && val == 0.0)
1153  count++;
1154  else if (type == L_GREATER_THAN_ZERO && val > 0.0)
1155  count++;
1156  }
1157 
1158  *pcount = count;
1159  return 0;
1160 }
1161 
1162 
1181 NUMA *
1183  l_int32 first,
1184  l_int32 last)
1185 {
1186 l_int32 n, i;
1187 l_float32 val, startx, delx;
1188 NUMA *nad;
1189 
1190  PROCNAME("numaClipToInterval");
1191 
1192  if (!nas)
1193  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
1194  if ((n = numaGetCount(nas)) == 0)
1195  return (NUMA *)ERROR_PTR("nas is empty", procName, NULL);
1196  if (first < 0 || first > last)
1197  return (NUMA *)ERROR_PTR("range not valid", procName, NULL);
1198  if (first >= n)
1199  return (NUMA *)ERROR_PTR("no elements in range", procName, NULL);
1200 
1201  last = L_MIN(last, n - 1);
1202  if ((nad = numaCreate(last - first + 1)) == NULL)
1203  return (NUMA *)ERROR_PTR("nad not made", procName, NULL);
1204  for (i = first; i <= last; i++) {
1205  numaGetFValue(nas, i, &val);
1206  numaAddNumber(nad, val);
1207  }
1208  numaGetParameters(nas, &startx, &delx);
1209  numaSetParameters(nad, startx + first * delx, delx);
1210  return nad;
1211 }
1212 
1213 
1230 NUMA *
1232  l_float32 thresh,
1233  l_int32 type)
1234 {
1235 l_int32 n, i, ival;
1236 l_float32 fval;
1237 NUMA *nai;
1238 
1239  PROCNAME("numaMakeThresholdIndicator");
1240 
1241  if (!nas)
1242  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
1243  if ((n = numaGetCount(nas)) == 0)
1244  return (NUMA *)ERROR_PTR("nas is empty", procName, NULL);
1245 
1246  nai = numaCreate(n);
1247  for (i = 0; i < n; i++) {
1248  numaGetFValue(nas, i, &fval);
1249  ival = 0;
1250  switch (type)
1251  {
1252  case L_SELECT_IF_LT:
1253  if (fval < thresh) ival = 1;
1254  break;
1255  case L_SELECT_IF_GT:
1256  if (fval > thresh) ival = 1;
1257  break;
1258  case L_SELECT_IF_LTE:
1259  if (fval <= thresh) ival = 1;
1260  break;
1261  case L_SELECT_IF_GTE:
1262  if (fval >= thresh) ival = 1;
1263  break;
1264  default:
1265  numaDestroy(&nai);
1266  return (NUMA *)ERROR_PTR("invalid type", procName, NULL);
1267  }
1268  numaAddNumber(nai, ival);
1269  }
1270 
1271  return nai;
1272 }
1273 
1274 
1288 NUMA *
1290  l_int32 nsamp)
1291 {
1292 l_int32 n, i, j, ileft, iright;
1293 l_float32 left, right, binsize, lfract, rfract, sum, startx, delx;
1294 l_float32 *array;
1295 NUMA *nad;
1296 
1297  PROCNAME("numaUniformSampling");
1298 
1299  if (!nas)
1300  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
1301  if ((n = numaGetCount(nas)) == 0)
1302  return (NUMA *)ERROR_PTR("nas is empty", procName, NULL);
1303  if (nsamp <= 0)
1304  return (NUMA *)ERROR_PTR("nsamp must be > 0", procName, NULL);
1305 
1306  nad = numaCreate(nsamp);
1307  array = numaGetFArray(nas, L_NOCOPY);
1308  binsize = (l_float32)n / (l_float32)nsamp;
1309  numaGetParameters(nas, &startx, &delx);
1310  numaSetParameters(nad, startx, binsize * delx);
1311  left = 0.0;
1312  for (i = 0; i < nsamp; i++) {
1313  sum = 0.0;
1314  right = left + binsize;
1315  ileft = (l_int32)left;
1316  lfract = 1.0 - left + ileft;
1317  if (lfract >= 1.0) /* on left bin boundary */
1318  lfract = 0.0;
1319  iright = (l_int32)right;
1320  rfract = right - iright;
1321  iright = L_MIN(iright, n - 1);
1322  if (ileft == iright) { /* both are within the same original sample */
1323  sum += (lfract + rfract - 1.0) * array[ileft];
1324  } else {
1325  if (lfract > 0.0001) /* left fraction */
1326  sum += lfract * array[ileft];
1327  if (rfract > 0.0001) /* right fraction */
1328  sum += rfract * array[iright];
1329  for (j = ileft + 1; j < iright; j++) /* entire pixels */
1330  sum += array[j];
1331  }
1332 
1333  numaAddNumber(nad, sum);
1334  left = right;
1335  }
1336  return nad;
1337 }
1338 
1339 
1354 NUMA *
1356  NUMA *nas)
1357 {
1358 l_int32 n, i;
1359 l_float32 val1, val2;
1360 
1361  PROCNAME("numaReverse");
1362 
1363  if (!nas)
1364  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
1365  if (nad && nas != nad)
1366  return (NUMA *)ERROR_PTR("nad defined but != nas", procName, NULL);
1367 
1368  n = numaGetCount(nas);
1369  if (nad) { /* in-place */
1370  for (i = 0; i < n / 2; i++) {
1371  numaGetFValue(nad, i, &val1);
1372  numaGetFValue(nad, n - i - 1, &val2);
1373  numaSetValue(nad, i, val2);
1374  numaSetValue(nad, n - i - 1, val1);
1375  }
1376  } else {
1377  nad = numaCreate(n);
1378  for (i = n - 1; i >= 0; i--) {
1379  numaGetFValue(nas, i, &val1);
1380  numaAddNumber(nad, val1);
1381  }
1382  }
1383 
1384  /* Reverse the startx and delx fields */
1385  nad->startx = nas->startx + (n - 1) * nas->delx;
1386  nad->delx = -nas->delx;
1387  return nad;
1388 }
1389 
1390 
1391 /*----------------------------------------------------------------------*
1392  * Signal feature extraction *
1393  *----------------------------------------------------------------------*/
1409 NUMA *
1411  l_float32 thresh,
1412  l_float32 maxn)
1413 {
1414 l_int32 n, i, inrun;
1415 l_float32 maxval, threshval, fval, startx, delx, x0, x1;
1416 NUMA *nad;
1417 
1418  PROCNAME("numaLowPassIntervals");
1419 
1420  if (!nas)
1421  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
1422  if ((n = numaGetCount(nas)) == 0)
1423  return (NUMA *)ERROR_PTR("nas is empty", procName, NULL);
1424  if (thresh < 0.0 || thresh > 1.0)
1425  return (NUMA *)ERROR_PTR("invalid thresh", procName, NULL);
1426 
1427  /* The input threshold is a fraction of the max.
1428  * The first entry in nad is the value of the max. */
1429  if (maxn == 0.0)
1430  numaGetMax(nas, &maxval, NULL);
1431  else
1432  maxval = maxn;
1433  numaGetParameters(nas, &startx, &delx);
1434  threshval = thresh * maxval;
1435  nad = numaCreate(0);
1436  numaAddNumber(nad, maxval);
1437 
1438  /* Write pairs of pts (x0, x1) for the intervals */
1439  inrun = FALSE;
1440  for (i = 0; i < n; i++) {
1441  numaGetFValue(nas, i, &fval);
1442  if (fval < threshval && inrun == FALSE) { /* start a new run */
1443  inrun = TRUE;
1444  x0 = startx + i * delx;
1445  } else if (fval > threshval && inrun == TRUE) { /* end the run */
1446  inrun = FALSE;
1447  x1 = startx + i * delx;
1448  numaAddNumber(nad, x0);
1449  numaAddNumber(nad, x1);
1450  }
1451  }
1452  if (inrun == TRUE) { /* must end the last run */
1453  x1 = startx + (n - 1) * delx;
1454  numaAddNumber(nad, x0);
1455  numaAddNumber(nad, x1);
1456  }
1457 
1458  return nad;
1459 }
1460 
1461 
1486 NUMA *
1488  l_float32 thresh1,
1489  l_float32 thresh2,
1490  l_float32 maxn)
1491 {
1492 l_int32 n, i, istart, inband, output, sign;
1493 l_int32 startbelow, below, above, belowlast, abovelast;
1494 l_float32 maxval, threshval1, threshval2, fval, startx, delx, x0, x1;
1495 NUMA *nad;
1496 
1497  PROCNAME("numaThresholdEdges");
1498 
1499  if (!nas)
1500  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
1501  if ((n = numaGetCount(nas)) == 0)
1502  return (NUMA *)ERROR_PTR("nas is empty", procName, NULL);
1503  if (thresh1 < 0.0 || thresh1 > 1.0 || thresh2 < 0.0 || thresh2 > 1.0)
1504  return (NUMA *)ERROR_PTR("invalid thresholds", procName, NULL);
1505  if (thresh2 < thresh1)
1506  return (NUMA *)ERROR_PTR("thresh2 < thresh1", procName, NULL);
1507 
1508  /* The input thresholds are fractions of the max.
1509  * The first entry in nad is the value of the max used
1510  * here for normalization. */
1511  if (maxn == 0.0)
1512  numaGetMax(nas, &maxval, NULL);
1513  else
1514  maxval = maxn;
1515  numaGetMax(nas, &maxval, NULL);
1516  numaGetParameters(nas, &startx, &delx);
1517  threshval1 = thresh1 * maxval;
1518  threshval2 = thresh2 * maxval;
1519  nad = numaCreate(0);
1520  numaAddNumber(nad, maxval);
1521 
1522  /* Write triplets of pts (x0, x1, sign) for the edges.
1523  * First make sure we start search from outside the band.
1524  * Only one of {belowlast, abovelast} is true. */
1525  for (i = 0; i < n; i++) {
1526  istart = i;
1527  numaGetFValue(nas, i, &fval);
1528  belowlast = (fval < threshval1) ? TRUE : FALSE;
1529  abovelast = (fval > threshval2) ? TRUE : FALSE;
1530  if (belowlast == TRUE || abovelast == TRUE)
1531  break;
1532  }
1533  if (istart == n) /* no intervals found */
1534  return nad;
1535 
1536  /* x0 and x1 can only be set from outside the edge.
1537  * They are the values just before entering the band,
1538  * and just after entering the band. We can jump through
1539  * the band, in which case they differ by one index in nas. */
1540  inband = FALSE;
1541  startbelow = belowlast; /* one of these is true */
1542  output = FALSE;
1543  x0 = startx + istart * delx;
1544  for (i = istart + 1; i < n; i++) {
1545  numaGetFValue(nas, i, &fval);
1546  below = (fval < threshval1) ? TRUE : FALSE;
1547  above = (fval > threshval2) ? TRUE : FALSE;
1548  if (!inband && belowlast && above) { /* full jump up */
1549  x1 = startx + i * delx;
1550  sign = 1;
1551  startbelow = FALSE; /* for the next transition */
1552  output = TRUE;
1553  } else if (!inband && abovelast && below) { /* full jump down */
1554  x1 = startx + i * delx;
1555  sign = -1;
1556  startbelow = TRUE; /* for the next transition */
1557  output = TRUE;
1558  } else if (inband && startbelow && above) { /* exit rising; success */
1559  x1 = startx + i * delx;
1560  sign = 1;
1561  inband = FALSE;
1562  startbelow = FALSE; /* for the next transition */
1563  output = TRUE;
1564  } else if (inband && !startbelow && below) {
1565  /* exit falling; success */
1566  x1 = startx + i * delx;
1567  sign = -1;
1568  inband = FALSE;
1569  startbelow = TRUE; /* for the next transition */
1570  output = TRUE;
1571  } else if (inband && !startbelow && above) { /* exit rising; failure */
1572  x0 = startx + i * delx;
1573  inband = FALSE;
1574  } else if (inband && startbelow && below) { /* exit falling; failure */
1575  x0 = startx + i * delx;
1576  inband = FALSE;
1577  } else if (!inband && !above && !below) { /* enter */
1578  inband = TRUE;
1579  startbelow = belowlast;
1580  } else if (!inband && (above || below)) { /* outside and remaining */
1581  x0 = startx + i * delx; /* update position */
1582  }
1583  belowlast = below;
1584  abovelast = above;
1585  if (output) { /* we have exited; save new x0 */
1586  numaAddNumber(nad, x0);
1587  numaAddNumber(nad, x1);
1588  numaAddNumber(nad, sign);
1589  output = FALSE;
1590  x0 = startx + i * delx;
1591  }
1592  }
1593 
1594  return nad;
1595 }
1596 
1597 
1607 l_int32
1609  l_int32 span,
1610  l_int32 *pstart,
1611  l_int32 *pend)
1612 {
1613 l_int32 n, nspans;
1614 
1615  PROCNAME("numaGetSpanValues");
1616 
1617  if (!na)
1618  return ERROR_INT("na not defined", procName, 1);
1619  if ((n = numaGetCount(na)) == 0)
1620  return ERROR_INT("na is empty", procName, 1);
1621  if (n % 2 != 1)
1622  return ERROR_INT("n is not odd", procName, 1);
1623  nspans = n / 2;
1624  if (nspans < 0 || span >= nspans)
1625  return ERROR_INT("invalid span", procName, 1);
1626 
1627  if (pstart) numaGetIValue(na, 2 * span + 1, pstart);
1628  if (pend) numaGetIValue(na, 2 * span + 2, pend);
1629  return 0;
1630 }
1631 
1632 
1644 l_int32
1646  l_int32 edge,
1647  l_int32 *pstart,
1648  l_int32 *pend,
1649  l_int32 *psign)
1650 {
1651 l_int32 n, nedges;
1652 
1653  PROCNAME("numaGetEdgeValues");
1654 
1655  if (!na)
1656  return ERROR_INT("na not defined", procName, 1);
1657  if ((n = numaGetCount(na)) == 0)
1658  return ERROR_INT("na is empty", procName, 1);
1659  if (n % 3 != 1)
1660  return ERROR_INT("n % 3 is not 1", procName, 1);
1661  nedges = (n - 1) / 3;
1662  if (edge < 0 || edge >= nedges)
1663  return ERROR_INT("invalid edge", procName, 1);
1664 
1665  if (pstart) numaGetIValue(na, 3 * edge + 1, pstart);
1666  if (pend) numaGetIValue(na, 3 * edge + 2, pend);
1667  if (psign) numaGetIValue(na, 3 * edge + 3, psign);
1668  return 0;
1669 }
1670 
1671 
1672 /*----------------------------------------------------------------------*
1673  * Interpolation *
1674  *----------------------------------------------------------------------*/
1702 l_ok
1703 numaInterpolateEqxVal(l_float32 startx,
1704  l_float32 deltax,
1705  NUMA *nay,
1706  l_int32 type,
1707  l_float32 xval,
1708  l_float32 *pyval)
1709 {
1710 l_int32 i, n, i1, i2, i3;
1711 l_float32 x1, x2, x3, fy1, fy2, fy3, d1, d2, d3, del, fi, maxx;
1712 l_float32 *fa;
1713 
1714  PROCNAME("numaInterpolateEqxVal");
1715 
1716  if (!pyval)
1717  return ERROR_INT("&yval not defined", procName, 1);
1718  *pyval = 0.0;
1719  if (!nay)
1720  return ERROR_INT("nay not defined", procName, 1);
1721  if (deltax <= 0.0)
1722  return ERROR_INT("deltax not > 0", procName, 1);
1723  if (type != L_LINEAR_INTERP && type != L_QUADRATIC_INTERP)
1724  return ERROR_INT("invalid interp type", procName, 1);
1725  if ((n = numaGetCount(nay)) < 2)
1726  return ERROR_INT("not enough points", procName, 1);
1727  if (type == L_QUADRATIC_INTERP && n == 2) {
1728  type = L_LINEAR_INTERP;
1729  L_WARNING("only 2 points; using linear interp\n", procName);
1730  }
1731  maxx = startx + deltax * (n - 1);
1732  if (xval < startx || xval > maxx)
1733  return ERROR_INT("xval is out of bounds", procName, 1);
1734 
1735  fa = numaGetFArray(nay, L_NOCOPY);
1736  fi = (xval - startx) / deltax;
1737  i = (l_int32)fi;
1738  del = fi - i;
1739  if (del == 0.0) { /* no interpolation required */
1740  *pyval = fa[i];
1741  return 0;
1742  }
1743 
1744  if (type == L_LINEAR_INTERP) {
1745  *pyval = fa[i] + del * (fa[i + 1] - fa[i]);
1746  return 0;
1747  }
1748 
1749  /* Quadratic interpolation */
1750  d1 = d3 = 0.5 / (deltax * deltax);
1751  d2 = -2. * d1;
1752  if (i == 0) {
1753  i1 = i;
1754  i2 = i + 1;
1755  i3 = i + 2;
1756  } else {
1757  i1 = i - 1;
1758  i2 = i;
1759  i3 = i + 1;
1760  }
1761  x1 = startx + i1 * deltax;
1762  x2 = startx + i2 * deltax;
1763  x3 = startx + i3 * deltax;
1764  fy1 = d1 * fa[i1];
1765  fy2 = d2 * fa[i2];
1766  fy3 = d3 * fa[i3];
1767  *pyval = fy1 * (xval - x2) * (xval - x3) +
1768  fy2 * (xval - x1) * (xval - x3) +
1769  fy3 * (xval - x1) * (xval - x2);
1770  return 0;
1771 }
1772 
1773 
1794 l_ok
1796  NUMA *nay,
1797  l_int32 type,
1798  l_float32 xval,
1799  l_float32 *pyval)
1800 {
1801 l_int32 i, im, nx, ny, i1, i2, i3;
1802 l_float32 delu, dell, fract, d1, d2, d3;
1803 l_float32 minx, maxx;
1804 l_float32 *fax, *fay;
1805 
1806  PROCNAME("numaInterpolateArbxVal");
1807 
1808  if (!pyval)
1809  return ERROR_INT("&yval not defined", procName, 1);
1810  *pyval = 0.0;
1811  if (!nax)
1812  return ERROR_INT("nax not defined", procName, 1);
1813  if (!nay)
1814  return ERROR_INT("nay not defined", procName, 1);
1815  if (type != L_LINEAR_INTERP && type != L_QUADRATIC_INTERP)
1816  return ERROR_INT("invalid interp type", procName, 1);
1817  ny = numaGetCount(nay);
1818  nx = numaGetCount(nax);
1819  if (nx != ny)
1820  return ERROR_INT("nax and nay not same size arrays", procName, 1);
1821  if (ny < 2)
1822  return ERROR_INT("not enough points", procName, 1);
1823  if (type == L_QUADRATIC_INTERP && ny == 2) {
1824  type = L_LINEAR_INTERP;
1825  L_WARNING("only 2 points; using linear interp\n", procName);
1826  }
1827  numaGetFValue(nax, 0, &minx);
1828  numaGetFValue(nax, nx - 1, &maxx);
1829  if (xval < minx || xval > maxx)
1830  return ERROR_INT("xval is out of bounds", procName, 1);
1831 
1832  fax = numaGetFArray(nax, L_NOCOPY);
1833  fay = numaGetFArray(nay, L_NOCOPY);
1834 
1835  /* Linear search for interval. We are guaranteed
1836  * to either return or break out of the loop.
1837  * In addition, we are assured that fax[i] - fax[im] > 0.0 */
1838  if (xval == fax[0]) {
1839  *pyval = fay[0];
1840  return 0;
1841  }
1842  im = 0;
1843  dell = 0.0;
1844  for (i = 1; i < nx; i++) {
1845  delu = fax[i] - xval;
1846  if (delu >= 0.0) { /* we've passed it */
1847  if (delu == 0.0) {
1848  *pyval = fay[i];
1849  return 0;
1850  }
1851  im = i - 1;
1852  dell = xval - fax[im]; /* >= 0 */
1853  break;
1854  }
1855  }
1856  fract = dell / (fax[i] - fax[im]);
1857 
1858  if (type == L_LINEAR_INTERP) {
1859  *pyval = fay[i] + fract * (fay[i + 1] - fay[i]);
1860  return 0;
1861  }
1862 
1863  /* Quadratic interpolation */
1864  if (im == 0) {
1865  i1 = im;
1866  i2 = im + 1;
1867  i3 = im + 2;
1868  } else {
1869  i1 = im - 1;
1870  i2 = im;
1871  i3 = im + 1;
1872  }
1873  d1 = (fax[i1] - fax[i2]) * (fax[i1] - fax[i3]);
1874  d2 = (fax[i2] - fax[i1]) * (fax[i2] - fax[i3]);
1875  d3 = (fax[i3] - fax[i1]) * (fax[i3] - fax[i2]);
1876  *pyval = fay[i1] * (xval - fax[i2]) * (xval - fax[i3]) / d1 +
1877  fay[i2] * (xval - fax[i1]) * (xval - fax[i3]) / d2 +
1878  fay[i3] * (xval - fax[i1]) * (xval - fax[i2]) / d3;
1879  return 0;
1880 }
1881 
1882 
1911 l_ok
1913  l_float32 deltax,
1914  NUMA *nasy,
1915  l_int32 type,
1916  l_float32 x0,
1917  l_float32 x1,
1918  l_int32 npts,
1919  NUMA **pnax,
1920  NUMA **pnay)
1921 {
1922 l_int32 i, n;
1923 l_float32 x, yval, maxx, delx;
1924 NUMA *nax, *nay;
1925 
1926  PROCNAME("numaInterpolateEqxInterval");
1927 
1928  if (pnax) *pnax = NULL;
1929  if (!pnay)
1930  return ERROR_INT("&nay not defined", procName, 1);
1931  *pnay = NULL;
1932  if (!nasy)
1933  return ERROR_INT("nasy not defined", procName, 1);
1934  if ((n = numaGetCount(nasy)) < 2)
1935  return ERROR_INT("n < 2", procName, 1);
1936  if (deltax <= 0.0)
1937  return ERROR_INT("deltax not > 0", procName, 1);
1938  if (type != L_LINEAR_INTERP && type != L_QUADRATIC_INTERP)
1939  return ERROR_INT("invalid interp type", procName, 1);
1940  if (type == L_QUADRATIC_INTERP && n == 2) {
1941  type = L_LINEAR_INTERP;
1942  L_WARNING("only 2 points; using linear interp\n", procName);
1943  }
1944  maxx = startx + deltax * (n - 1);
1945  if (x0 < startx || x1 > maxx || x1 <= x0)
1946  return ERROR_INT("[x0 ... x1] is not valid", procName, 1);
1947  if (npts < 3)
1948  return ERROR_INT("npts < 3", procName, 1);
1949  delx = (x1 - x0) / (l_float32)(npts - 1); /* delx is for output nay */
1950 
1951  if ((nay = numaCreate(npts)) == NULL)
1952  return ERROR_INT("nay not made", procName, 1);
1953  numaSetParameters(nay, x0, delx);
1954  *pnay = nay;
1955  if (pnax) {
1956  nax = numaCreate(npts);
1957  *pnax = nax;
1958  }
1959 
1960  for (i = 0; i < npts; i++) {
1961  x = x0 + i * delx;
1962  if (pnax)
1963  numaAddNumber(nax, x);
1964  numaInterpolateEqxVal(startx, deltax, nasy, type, x, &yval);
1965  numaAddNumber(nay, yval);
1966  }
1967 
1968  return 0;
1969 }
1970 
1971 
2000 l_ok
2002  NUMA *nay,
2003  l_int32 type,
2004  l_float32 x0,
2005  l_float32 x1,
2006  l_int32 npts,
2007  NUMA **pnadx,
2008  NUMA **pnady)
2009 {
2010 l_int32 i, im, j, nx, ny, i1, i2, i3, sorted;
2011 l_int32 *index;
2012 l_float32 del, xval, yval, excess, fract, minx, maxx, d1, d2, d3;
2013 l_float32 *fax, *fay;
2014 NUMA *nasx, *nasy, *nadx, *nady;
2015 
2016  PROCNAME("numaInterpolateArbxInterval");
2017 
2018  if (pnadx) *pnadx = NULL;
2019  if (!pnady)
2020  return ERROR_INT("&nady not defined", procName, 1);
2021  *pnady = NULL;
2022  if (!nay)
2023  return ERROR_INT("nay not defined", procName, 1);
2024  if (!nax)
2025  return ERROR_INT("nax not defined", procName, 1);
2026  if (type != L_LINEAR_INTERP && type != L_QUADRATIC_INTERP)
2027  return ERROR_INT("invalid interp type", procName, 1);
2028  if (x0 > x1)
2029  return ERROR_INT("x0 > x1", procName, 1);
2030  ny = numaGetCount(nay);
2031  nx = numaGetCount(nax);
2032  if (nx != ny)
2033  return ERROR_INT("nax and nay not same size arrays", procName, 1);
2034  if (ny < 2)
2035  return ERROR_INT("not enough points", procName, 1);
2036  if (type == L_QUADRATIC_INTERP && ny == 2) {
2037  type = L_LINEAR_INTERP;
2038  L_WARNING("only 2 points; using linear interp\n", procName);
2039  }
2040  numaGetMin(nax, &minx, NULL);
2041  numaGetMax(nax, &maxx, NULL);
2042  if (x0 < minx || x1 > maxx)
2043  return ERROR_INT("xval is out of bounds", procName, 1);
2044 
2045  /* Make sure that nax is sorted in increasing order */
2046  numaIsSorted(nax, L_SORT_INCREASING, &sorted);
2047  if (!sorted) {
2048  L_WARNING("we are sorting nax in increasing order\n", procName);
2049  numaSortPair(nax, nay, L_SORT_INCREASING, &nasx, &nasy);
2050  } else {
2051  nasx = numaClone(nax);
2052  nasy = numaClone(nay);
2053  }
2054 
2055  fax = numaGetFArray(nasx, L_NOCOPY);
2056  fay = numaGetFArray(nasy, L_NOCOPY);
2057 
2058  /* Get array of indices into fax for interpolated locations */
2059  if ((index = (l_int32 *)LEPT_CALLOC(npts, sizeof(l_int32))) == NULL) {
2060  numaDestroy(&nasx);
2061  numaDestroy(&nasy);
2062  return ERROR_INT("ind not made", procName, 1);
2063  }
2064  del = (x1 - x0) / (npts - 1.0);
2065  for (i = 0, j = 0; j < nx && i < npts; i++) {
2066  xval = x0 + i * del;
2067  while (j < nx - 1 && xval > fax[j])
2068  j++;
2069  if (xval == fax[j])
2070  index[i] = L_MIN(j, nx - 1);
2071  else /* the index of fax[] is just below xval */
2072  index[i] = L_MAX(j - 1, 0);
2073  }
2074 
2075  /* For each point to be interpolated, get the y value */
2076  nady = numaCreate(npts);
2077  *pnady = nady;
2078  if (pnadx) {
2079  nadx = numaCreate(npts);
2080  *pnadx = nadx;
2081  }
2082  for (i = 0; i < npts; i++) {
2083  xval = x0 + i * del;
2084  if (pnadx)
2085  numaAddNumber(nadx, xval);
2086  im = index[i];
2087  excess = xval - fax[im];
2088  if (excess == 0.0) {
2089  numaAddNumber(nady, fay[im]);
2090  continue;
2091  }
2092  fract = excess / (fax[im + 1] - fax[im]);
2093 
2094  if (type == L_LINEAR_INTERP) {
2095  yval = fay[im] + fract * (fay[im + 1] - fay[im]);
2096  numaAddNumber(nady, yval);
2097  continue;
2098  }
2099 
2100  /* Quadratic interpolation */
2101  if (im == 0) {
2102  i1 = im;
2103  i2 = im + 1;
2104  i3 = im + 2;
2105  } else {
2106  i1 = im - 1;
2107  i2 = im;
2108  i3 = im + 1;
2109  }
2110  d1 = (fax[i1] - fax[i2]) * (fax[i1] - fax[i3]);
2111  d2 = (fax[i2] - fax[i1]) * (fax[i2] - fax[i3]);
2112  d3 = (fax[i3] - fax[i1]) * (fax[i3] - fax[i2]);
2113  yval = fay[i1] * (xval - fax[i2]) * (xval - fax[i3]) / d1 +
2114  fay[i2] * (xval - fax[i1]) * (xval - fax[i3]) / d2 +
2115  fay[i3] * (xval - fax[i1]) * (xval - fax[i2]) / d3;
2116  numaAddNumber(nady, yval);
2117  }
2118 
2119  LEPT_FREE(index);
2120  numaDestroy(&nasx);
2121  numaDestroy(&nasy);
2122  return 0;
2123 }
2124 
2125 
2126 /*----------------------------------------------------------------------*
2127  * Functions requiring interpolation *
2128  *----------------------------------------------------------------------*/
2161 l_ok
2163  l_float32 *pmaxval,
2164  NUMA *naloc,
2165  l_float32 *pmaxloc)
2166 {
2167 l_float32 val;
2168 l_float32 smaxval; /* start value of maximum sample, before interpolating */
2169 l_int32 n, imaxloc;
2170 l_float32 x1, x2, x3, y1, y2, y3, c1, c2, c3, a, b, xmax, ymax;
2171 
2172  PROCNAME("numaFitMax");
2173 
2174  if (pmaxval) *pmaxval = 0.0;
2175  if (pmaxloc) *pmaxloc = 0.0;
2176  if (!na)
2177  return ERROR_INT("na not defined", procName, 1);
2178  if ((n = numaGetCount(na)) == 0)
2179  return ERROR_INT("na is empty", procName, 1);
2180  if (!pmaxval)
2181  return ERROR_INT("&maxval not defined", procName, 1);
2182  if (!pmaxloc)
2183  return ERROR_INT("&maxloc not defined", procName, 1);
2184  if (naloc) {
2185  if (n != numaGetCount(naloc))
2186  return ERROR_INT("na and naloc of unequal size", procName, 1);
2187  }
2188 
2189  numaGetMax(na, &smaxval, &imaxloc);
2190 
2191  /* Simple case: max is at end point */
2192  if (imaxloc == 0 || imaxloc == n - 1) {
2193  *pmaxval = smaxval;
2194  if (naloc) {
2195  numaGetFValue(naloc, imaxloc, &val);
2196  *pmaxloc = val;
2197  } else {
2198  *pmaxloc = imaxloc;
2199  }
2200  return 0;
2201  }
2202 
2203  /* Interior point; use quadratic interpolation */
2204  y2 = smaxval;
2205  numaGetFValue(na, imaxloc - 1, &val);
2206  y1 = val;
2207  numaGetFValue(na, imaxloc + 1, &val);
2208  y3 = val;
2209  if (naloc) {
2210  numaGetFValue(naloc, imaxloc - 1, &val);
2211  x1 = val;
2212  numaGetFValue(naloc, imaxloc, &val);
2213  x2 = val;
2214  numaGetFValue(naloc, imaxloc + 1, &val);
2215  x3 = val;
2216  } else {
2217  x1 = imaxloc - 1;
2218  x2 = imaxloc;
2219  x3 = imaxloc + 1;
2220  }
2221 
2222  /* Can't interpolate; just use the max val in na
2223  * and the corresponding one in naloc */
2224  if (x1 == x2 || x1 == x3 || x2 == x3) {
2225  *pmaxval = y2;
2226  *pmaxloc = x2;
2227  return 0;
2228  }
2229 
2230  /* Use lagrangian interpolation; set dy/dx = 0 */
2231  c1 = y1 / ((x1 - x2) * (x1 - x3));
2232  c2 = y2 / ((x2 - x1) * (x2 - x3));
2233  c3 = y3 / ((x3 - x1) * (x3 - x2));
2234  a = c1 + c2 + c3;
2235  b = c1 * (x2 + x3) + c2 * (x1 + x3) + c3 * (x1 + x2);
2236  xmax = b / (2 * a);
2237  ymax = c1 * (xmax - x2) * (xmax - x3) +
2238  c2 * (xmax - x1) * (xmax - x3) +
2239  c3 * (xmax - x1) * (xmax - x2);
2240  *pmaxval = ymax;
2241  *pmaxloc = xmax;
2242 
2243  return 0;
2244 }
2245 
2246 
2267 l_ok
2269  NUMA *nay,
2270  l_float32 x0,
2271  l_float32 x1,
2272  l_int32 npts,
2273  NUMA **pnadx,
2274  NUMA **pnady)
2275 {
2276 l_int32 i, nx, ny;
2277 l_float32 minx, maxx, der, invdel;
2278 l_float32 *fay;
2279 NUMA *nady, *naiy;
2280 
2281  PROCNAME("numaDifferentiateInterval");
2282 
2283  if (pnadx) *pnadx = NULL;
2284  if (!pnady)
2285  return ERROR_INT("&nady not defined", procName, 1);
2286  *pnady = NULL;
2287  if (!nay)
2288  return ERROR_INT("nay not defined", procName, 1);
2289  if (!nax)
2290  return ERROR_INT("nax not defined", procName, 1);
2291  if (x0 > x1)
2292  return ERROR_INT("x0 > x1", procName, 1);
2293  ny = numaGetCount(nay);
2294  nx = numaGetCount(nax);
2295  if (nx != ny)
2296  return ERROR_INT("nax and nay not same size arrays", procName, 1);
2297  if (ny < 2)
2298  return ERROR_INT("not enough points", procName, 1);
2299  numaGetMin(nax, &minx, NULL);
2300  numaGetMax(nax, &maxx, NULL);
2301  if (x0 < minx || x1 > maxx)
2302  return ERROR_INT("xval is out of bounds", procName, 1);
2303  if (npts < 2)
2304  return ERROR_INT("npts < 2", procName, 1);
2305 
2306  /* Generate interpolated array over specified interval */
2307  if (numaInterpolateArbxInterval(nax, nay, L_LINEAR_INTERP, x0, x1,
2308  npts, pnadx, &naiy))
2309  return ERROR_INT("interpolation failed", procName, 1);
2310 
2311  nady = numaCreate(npts);
2312  *pnady = nady;
2313  invdel = 0.5 * ((l_float32)npts - 1.0) / (x1 - x0);
2314  fay = numaGetFArray(naiy, L_NOCOPY);
2315 
2316  /* Compute and save derivatives */
2317  der = 0.5 * invdel * (fay[1] - fay[0]);
2318  numaAddNumber(nady, der);
2319  for (i = 1; i < npts - 1; i++) {
2320  der = invdel * (fay[i + 1] - fay[i - 1]);
2321  numaAddNumber(nady, der);
2322  }
2323  der = 0.5 * invdel * (fay[npts - 1] - fay[npts - 2]);
2324  numaAddNumber(nady, der);
2325 
2326  numaDestroy(&naiy);
2327  return 0;
2328 }
2329 
2330 
2350 l_ok
2352  NUMA *nay,
2353  l_float32 x0,
2354  l_float32 x1,
2355  l_int32 npts,
2356  l_float32 *psum)
2357 {
2358 l_int32 i, nx, ny;
2359 l_float32 minx, maxx, sum, del;
2360 l_float32 *fay;
2361 NUMA *naiy;
2362 
2363  PROCNAME("numaIntegrateInterval");
2364 
2365  if (!psum)
2366  return ERROR_INT("&sum not defined", procName, 1);
2367  *psum = 0.0;
2368  if (!nay)
2369  return ERROR_INT("nay not defined", procName, 1);
2370  if (!nax)
2371  return ERROR_INT("nax not defined", procName, 1);
2372  if (x0 > x1)
2373  return ERROR_INT("x0 > x1", procName, 1);
2374  if (npts < 2)
2375  return ERROR_INT("npts < 2", procName, 1);
2376  ny = numaGetCount(nay);
2377  nx = numaGetCount(nax);
2378  if (nx != ny)
2379  return ERROR_INT("nax and nay not same size arrays", procName, 1);
2380  if (ny < 2)
2381  return ERROR_INT("not enough points", procName, 1);
2382  numaGetMin(nax, &minx, NULL);
2383  numaGetMax(nax, &maxx, NULL);
2384  if (x0 < minx || x1 > maxx)
2385  return ERROR_INT("xval is out of bounds", procName, 1);
2386 
2387  /* Generate interpolated array over specified interval */
2388  if (numaInterpolateArbxInterval(nax, nay, L_LINEAR_INTERP, x0, x1,
2389  npts, NULL, &naiy))
2390  return ERROR_INT("interpolation failed", procName, 1);
2391 
2392  del = (x1 - x0) / ((l_float32)npts - 1.0);
2393  fay = numaGetFArray(naiy, L_NOCOPY);
2394 
2395  /* Compute integral (simple trapezoid) */
2396  sum = 0.5 * (fay[0] + fay[npts - 1]);
2397  for (i = 1; i < npts - 1; i++)
2398  sum += fay[i];
2399  *psum = del * sum;
2400 
2401  numaDestroy(&naiy);
2402  return 0;
2403 }
2404 
2405 
2406 /*----------------------------------------------------------------------*
2407  * Sorting *
2408  *----------------------------------------------------------------------*/
2455 l_ok
2457  NUMA **pnasort,
2458  NUMA **pnaindex,
2459  NUMA **pnainvert,
2460  l_int32 sortorder,
2461  l_int32 sorttype)
2462 {
2463 l_int32 isize;
2464 l_float32 size;
2465 NUMA *naindex = NULL;
2466 
2467  PROCNAME("numaSortGeneral");
2468 
2469  if (pnasort) *pnasort = NULL;
2470  if (pnaindex) *pnaindex = NULL;
2471  if (pnainvert) *pnainvert = NULL;
2472  if (!na)
2473  return ERROR_INT("na not defined", procName, 1);
2474  if (sortorder != L_SORT_INCREASING && sortorder != L_SORT_DECREASING)
2475  return ERROR_INT("invalid sort order", procName, 1);
2476  if (sorttype != L_SHELL_SORT && sorttype != L_BIN_SORT)
2477  return ERROR_INT("invalid sort type", procName, 1);
2478  if (!pnasort && !pnaindex && !pnainvert)
2479  return ERROR_INT("nothing to do", procName, 1);
2480 
2481  if (sorttype == L_BIN_SORT) {
2482  numaGetMax(na, &size, NULL);
2483  isize = (l_int32)size;
2484  if (isize > MaxInitPtraSize - 1) {
2485  L_WARNING("array too large; using shell sort\n", procName);
2486  sorttype = L_SHELL_SORT;
2487  } else {
2488  naindex = numaGetBinSortIndex(na, sortorder);
2489  }
2490  }
2491 
2492  if (sorttype == L_SHELL_SORT)
2493  naindex = numaGetSortIndex(na, sortorder);
2494 
2495  if (pnasort)
2496  *pnasort = numaSortByIndex(na, naindex);
2497  if (pnainvert)
2498  *pnainvert = numaInvertMap(naindex);
2499  if (pnaindex)
2500  *pnaindex = naindex;
2501  else
2502  numaDestroy(&naindex);
2503  return 0;
2504 }
2505 
2506 
2520 NUMA *
2522  l_int32 sortorder)
2523 {
2524 l_int32 type;
2525 
2526  PROCNAME("numaSortAutoSelect");
2527 
2528  if (!nas)
2529  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
2530  if (numaGetCount(nas) == 0) {
2531  L_WARNING("nas is empty; returning copy\n", procName);
2532  return numaCopy(nas);
2533  }
2534  if (sortorder != L_SORT_INCREASING && sortorder != L_SORT_DECREASING)
2535  return (NUMA *)ERROR_PTR("invalid sort order", procName, NULL);
2536 
2537  type = numaChooseSortType(nas);
2538  if (type != L_SHELL_SORT && type != L_BIN_SORT)
2539  return (NUMA *)ERROR_PTR("invalid sort type", procName, NULL);
2540 
2541  if (type == L_BIN_SORT)
2542  return numaBinSort(nas, sortorder);
2543  else /* shell sort */
2544  return numaSort(NULL, nas, sortorder);
2545 }
2546 
2547 
2561 NUMA *
2563  l_int32 sortorder)
2564 {
2565 l_int32 type;
2566 
2567  PROCNAME("numaSortIndexAutoSelect");
2568 
2569  if (!nas)
2570  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
2571  if (numaGetCount(nas) == 0) {
2572  L_WARNING("nas is empty; returning copy\n", procName);
2573  return numaCopy(nas);
2574  }
2575  if (sortorder != L_SORT_INCREASING && sortorder != L_SORT_DECREASING)
2576  return (NUMA *)ERROR_PTR("invalid sort order", procName, NULL);
2577  type = numaChooseSortType(nas);
2578  if (type != L_SHELL_SORT && type != L_BIN_SORT)
2579  return (NUMA *)ERROR_PTR("invalid sort type", procName, NULL);
2580 
2581  if (type == L_BIN_SORT)
2582  return numaGetBinSortIndex(nas, sortorder);
2583  else /* shell sort */
2584  return numaGetSortIndex(nas, sortorder);
2585 }
2586 
2587 
2601 l_int32
2603 {
2604 l_int32 n;
2605 l_float32 minval, maxval;
2606 
2607  PROCNAME("numaChooseSortType");
2608 
2609  if (!nas)
2610  return ERROR_INT("nas not defined", procName, UNDEF);
2611 
2612  /* If small histogram or negative values; use shell sort */
2613  numaGetMin(nas, &minval, NULL);
2614  n = numaGetCount(nas);
2615  if (minval < 0.0 || n < 200)
2616  return L_SHELL_SORT;
2617 
2618  /* If large maxval, use shell sort */
2619  numaGetMax(nas, &maxval, NULL);
2620  if (maxval > MaxInitPtraSize - 1)
2621  return L_SHELL_SORT;
2622 
2623  /* Otherwise, need to compare nlog(n) with maxval.
2624  * The factor of 0.003 was determined by comparing times for
2625  * different histogram sizes and maxval. It is very small
2626  * because binsort is fast and shell sort gets slow for large n. */
2627  if (n * log((l_float32)n) < 0.003 * maxval)
2628  return L_SHELL_SORT;
2629  else
2630  return L_BIN_SORT;
2631 }
2632 
2633 
2649 NUMA *
2651  NUMA *nain,
2652  l_int32 sortorder)
2653 {
2654 l_int32 i, n, gap, j;
2655 l_float32 tmp;
2656 l_float32 *array;
2657 
2658  PROCNAME("numaSort");
2659 
2660  if (!nain)
2661  return (NUMA *)ERROR_PTR("nain not defined", procName, NULL);
2662  if (sortorder != L_SORT_INCREASING && sortorder != L_SORT_DECREASING)
2663  return (NUMA *)ERROR_PTR("invalid sort order", procName, NULL);
2664 
2665  /* Make naout if necessary; otherwise do in-place */
2666  if (!naout)
2667  naout = numaCopy(nain);
2668  else if (nain != naout)
2669  return (NUMA *)ERROR_PTR("invalid: not in-place", procName, NULL);
2670  if ((n = numaGetCount(naout)) == 0) {
2671  L_WARNING("naout is empty\n", procName);
2672  return naout;
2673  }
2674  array = naout->array; /* operate directly on the array */
2675  n = numaGetCount(naout);
2676 
2677  /* Shell sort */
2678  for (gap = n/2; gap > 0; gap = gap / 2) {
2679  for (i = gap; i < n; i++) {
2680  for (j = i - gap; j >= 0; j -= gap) {
2681  if ((sortorder == L_SORT_INCREASING &&
2682  array[j] > array[j + gap]) ||
2683  (sortorder == L_SORT_DECREASING &&
2684  array[j] < array[j + gap]))
2685  {
2686  tmp = array[j];
2687  array[j] = array[j + gap];
2688  array[j + gap] = tmp;
2689  }
2690  }
2691  }
2692  }
2693 
2694  return naout;
2695 }
2696 
2697 
2717 NUMA *
2719  l_int32 sortorder)
2720 {
2721 NUMA *nat, *nad;
2722 
2723  PROCNAME("numaBinSort");
2724 
2725  if (!nas)
2726  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
2727  if (numaGetCount(nas) == 0) {
2728  L_WARNING("nas is empty; returning copy\n", procName);
2729  return numaCopy(nas);
2730  }
2731  if (sortorder != L_SORT_INCREASING && sortorder != L_SORT_DECREASING)
2732  return (NUMA *)ERROR_PTR("invalid sort order", procName, NULL);
2733 
2734  if ((nat = numaGetBinSortIndex(nas, sortorder)) == NULL)
2735  return (NUMA *)ERROR_PTR("bin sort failed", procName, NULL);
2736  nad = numaSortByIndex(nas, nat);
2737  numaDestroy(&nat);
2738  return nad;
2739 }
2740 
2741 
2750 NUMA *
2752  l_int32 sortorder)
2753 {
2754 l_int32 i, n, gap, j;
2755 l_float32 tmp;
2756 l_float32 *array; /* copy of input array */
2757 l_float32 *iarray; /* array of indices */
2758 NUMA *naisort;
2759 
2760  PROCNAME("numaGetSortIndex");
2761 
2762  if (!na)
2763  return (NUMA *)ERROR_PTR("na not defined", procName, NULL);
2764  if (numaGetCount(na) == 0) {
2765  L_WARNING("na is empty\n", procName);
2766  return numaCreate(1);
2767  }
2768  if (sortorder != L_SORT_INCREASING && sortorder != L_SORT_DECREASING)
2769  return (NUMA *)ERROR_PTR("invalid sortorder", procName, NULL);
2770 
2771  n = numaGetCount(na);
2772  if ((array = numaGetFArray(na, L_COPY)) == NULL)
2773  return (NUMA *)ERROR_PTR("array not made", procName, NULL);
2774  if ((iarray = (l_float32 *)LEPT_CALLOC(n, sizeof(l_float32))) == NULL) {
2775  LEPT_FREE(array);
2776  return (NUMA *)ERROR_PTR("iarray not made", procName, NULL);
2777  }
2778  for (i = 0; i < n; i++)
2779  iarray[i] = i;
2780 
2781  /* Shell sort */
2782  for (gap = n/2; gap > 0; gap = gap / 2) {
2783  for (i = gap; i < n; i++) {
2784  for (j = i - gap; j >= 0; j -= gap) {
2785  if ((sortorder == L_SORT_INCREASING &&
2786  array[j] > array[j + gap]) ||
2787  (sortorder == L_SORT_DECREASING &&
2788  array[j] < array[j + gap]))
2789  {
2790  tmp = array[j];
2791  array[j] = array[j + gap];
2792  array[j + gap] = tmp;
2793  tmp = iarray[j];
2794  iarray[j] = iarray[j + gap];
2795  iarray[j + gap] = tmp;
2796  }
2797  }
2798  }
2799  }
2800 
2801  naisort = numaCreate(n);
2802  for (i = 0; i < n; i++)
2803  numaAddNumber(naisort, iarray[i]);
2804 
2805  LEPT_FREE(array);
2806  LEPT_FREE(iarray);
2807  return naisort;
2808 }
2809 
2810 
2832 NUMA *
2834  l_int32 sortorder)
2835 {
2836 l_int32 i, n, isize, ival, imax;
2837 l_float32 minsize, size;
2838 NUMA *na, *nai, *nad;
2839 L_PTRA *paindex;
2840 
2841  PROCNAME("numaGetBinSortIndex");
2842 
2843  if (!nas)
2844  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
2845  if (numaGetCount(nas) == 0) {
2846  L_WARNING("nas is empty\n", procName);
2847  return numaCreate(1);
2848  }
2849  if (sortorder != L_SORT_INCREASING && sortorder != L_SORT_DECREASING)
2850  return (NUMA *)ERROR_PTR("invalid sort order", procName, NULL);
2851  numaGetMin(nas, &minsize, NULL);
2852  if (minsize < 0)
2853  return (NUMA *)ERROR_PTR("nas has negative numbers", procName, NULL);
2854  numaGetMax(nas, &size, NULL);
2855  isize = (l_int32)size;
2856  if (isize > MaxInitPtraSize - 1) {
2857  L_ERROR("array too large: %d elements > max size = %d\n",
2858  procName, isize, MaxInitPtraSize - 1);
2859  return NULL;
2860  }
2861 
2862  /* Set up a ptra holding numa at indices for which there
2863  * are values in nas. Suppose nas has the value 230 at index
2864  * 7355. A numa holding the index 7355 is created and stored
2865  * at the ptra index 230. If there is another value of 230
2866  * in nas, its index is added to the same numa (at index 230
2867  * in the ptra). When finished, the ptra can be scanned for numa,
2868  * and the original indices in the nas can be read out. In this
2869  * way, the ptra effectively sorts the input numbers in the nas. */
2870  paindex = ptraCreate(isize + 1);
2871  n = numaGetCount(nas);
2872  for (i = 0; i < n; i++) {
2873  numaGetIValue(nas, i, &ival);
2874  nai = (NUMA *)ptraGetPtrToItem(paindex, ival);
2875  if (!nai) { /* make it; no shifting will occur */
2876  nai = numaCreate(1);
2877  ptraInsert(paindex, ival, nai, L_MIN_DOWNSHIFT);
2878  }
2879  numaAddNumber(nai, i);
2880  }
2881 
2882  /* Sort by scanning the ptra, extracting numas and pulling
2883  * the (index into nas) numbers out of each numa, taken
2884  * successively in requested order. */
2885  ptraGetMaxIndex(paindex, &imax);
2886  nad = numaCreate(0);
2887  if (sortorder == L_SORT_INCREASING) {
2888  for (i = 0; i <= imax; i++) {
2889  na = (NUMA *)ptraRemove(paindex, i, L_NO_COMPACTION);
2890  if (!na) continue;
2891  numaJoin(nad, na, 0, -1);
2892  numaDestroy(&na);
2893  }
2894  } else { /* L_SORT_DECREASING */
2895  for (i = imax; i >= 0; i--) {
2896  na = (NUMA *)ptraRemoveLast(paindex);
2897  if (!na) break; /* they've all been removed */
2898  numaJoin(nad, na, 0, -1);
2899  numaDestroy(&na);
2900  }
2901  }
2902 
2903  ptraDestroy(&paindex, FALSE, FALSE);
2904  return nad;
2905 }
2906 
2907 
2915 NUMA *
2917  NUMA *naindex)
2918 {
2919 l_int32 i, n, ni, index;
2920 l_float32 val;
2921 NUMA *nad;
2922 
2923  PROCNAME("numaSortByIndex");
2924 
2925  if (!nas)
2926  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
2927  if (!naindex)
2928  return (NUMA *)ERROR_PTR("naindex not defined", procName, NULL);
2929  n = numaGetCount(nas);
2930  ni = numaGetCount(naindex);
2931  if (n != ni)
2932  return (NUMA *)ERROR_PTR("numa sizes differ", procName, NULL);
2933  if (n == 0) {
2934  L_WARNING("nas is empty\n", procName);
2935  return numaCopy(nas);
2936  }
2937 
2938  nad = numaCreate(n);
2939  for (i = 0; i < n; i++) {
2940  numaGetIValue(naindex, i, &index);
2941  numaGetFValue(nas, index, &val);
2942  numaAddNumber(nad, val);
2943  }
2944 
2945  return nad;
2946 }
2947 
2948 
2964 l_int32
2966  l_int32 sortorder,
2967  l_int32 *psorted)
2968 {
2969 l_int32 i, n;
2970 l_float32 prevval, val;
2971 
2972  PROCNAME("numaIsSorted");
2973 
2974  if (!psorted)
2975  return ERROR_INT("&sorted not defined", procName, 1);
2976  *psorted = FALSE;
2977  if (!nas)
2978  return ERROR_INT("nas not defined", procName, 1);
2979  if ((n = numaGetCount(nas))== 0) {
2980  L_WARNING("nas is empty\n", procName);
2981  *psorted = TRUE;
2982  return 0;
2983  }
2984  if (sortorder != L_SORT_INCREASING && sortorder != L_SORT_DECREASING)
2985  return ERROR_INT("invalid sortorder", procName, 1);
2986 
2987  n = numaGetCount(nas);
2988  numaGetFValue(nas, 0, &prevval);
2989  for (i = 1; i < n; i++) {
2990  numaGetFValue(nas, i, &val);
2991  if ((sortorder == L_SORT_INCREASING && val < prevval) ||
2992  (sortorder == L_SORT_DECREASING && val > prevval))
2993  return 0;
2994  }
2995 
2996  *psorted = TRUE;
2997  return 0;
2998 }
2999 
3000 
3016 l_ok
3018  NUMA *nay,
3019  l_int32 sortorder,
3020  NUMA **pnasx,
3021  NUMA **pnasy)
3022 {
3023 l_int32 sorted;
3024 NUMA *naindex;
3025 
3026  PROCNAME("numaSortPair");
3027 
3028  if (pnasx) *pnasx = NULL;
3029  if (pnasy) *pnasy = NULL;
3030  if (!pnasx || !pnasy)
3031  return ERROR_INT("&nasx and/or &nasy not defined", procName, 1);
3032  if (!nax)
3033  return ERROR_INT("nax not defined", procName, 1);
3034  if (!nay)
3035  return ERROR_INT("nay not defined", procName, 1);
3036  if (sortorder != L_SORT_INCREASING && sortorder != L_SORT_DECREASING)
3037  return ERROR_INT("invalid sortorder", procName, 1);
3038 
3039  numaIsSorted(nax, sortorder, &sorted);
3040  if (sorted == TRUE) {
3041  *pnasx = numaCopy(nax);
3042  *pnasy = numaCopy(nay);
3043  } else {
3044  naindex = numaGetSortIndex(nax, sortorder);
3045  *pnasx = numaSortByIndex(nax, naindex);
3046  *pnasy = numaSortByIndex(nay, naindex);
3047  numaDestroy(&naindex);
3048  }
3049 
3050  return 0;
3051 }
3052 
3053 
3067 NUMA *
3069 {
3070 l_int32 i, n, val, error;
3071 l_int32 *test;
3072 NUMA *nad;
3073 
3074  PROCNAME("numaInvertMap");
3075 
3076  if (!nas)
3077  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
3078  if ((n = numaGetCount(nas)) == 0) {
3079  L_WARNING("nas is empty\n", procName);
3080  return numaCopy(nas);
3081  }
3082 
3083  nad = numaMakeConstant(0.0, n);
3084  test = (l_int32 *)LEPT_CALLOC(n, sizeof(l_int32));
3085  error = 0;
3086  for (i = 0; i < n; i++) {
3087  numaGetIValue(nas, i, &val);
3088  if (val >= n) {
3089  error = 1;
3090  break;
3091  }
3092  numaReplaceNumber(nad, val, i);
3093  if (test[val] == 0) {
3094  test[val] = 1;
3095  } else {
3096  error = 1;
3097  break;
3098  }
3099  }
3100 
3101  LEPT_FREE(test);
3102  if (error) {
3103  numaDestroy(&nad);
3104  return (NUMA *)ERROR_PTR("nas not invertible", procName, NULL);
3105  }
3106 
3107  return nad;
3108 }
3109 
3123 l_ok
3125  l_float32 val)
3126 {
3127 l_int32 index;
3128 
3129  PROCNAME("numaAddSorted");
3130 
3131  if (!na)
3132  return ERROR_INT("na not defined", procName, 1);
3133 
3134  if (numaFindSortedLoc(na, val, &index) == 1)
3135  return ERROR_INT("insert failure", procName, 1);
3136  numaInsertNumber(na, index, val);
3137  return 0;
3138 }
3139 
3140 
3163 l_ok
3165  l_float32 val,
3166  l_int32 *pindex)
3167 {
3168 l_int32 n, increasing, lindex, rindex, midindex;
3169 l_float32 val0, valn, valmid;
3170 
3171  PROCNAME("numaFindSortedLoc");
3172 
3173  if (!pindex)
3174  return ERROR_INT("&index not defined", procName, 1);
3175  *pindex = 0;
3176  if (!na)
3177  return ERROR_INT("na not defined", procName, 1);
3178 
3179  n = numaGetCount(na);
3180  if (n == 0) return 0;
3181  numaGetFValue(na, 0, &val0);
3182  if (n == 1) { /* use increasing sort order */
3183  if (val >= val0)
3184  *pindex = 1;
3185  return 0;
3186  }
3187 
3188  /* ----------------- n >= 2 ----------------- */
3189  numaGetFValue(na, n - 1, &valn);
3190  increasing = (valn >= val0) ? 1 : 0; /* sort order */
3191 
3192  /* Check if outside bounds of existing array */
3193  if (increasing) {
3194  if (val < val0) {
3195  *pindex = 0;
3196  return 0;
3197  } else if (val > valn) {
3198  *pindex = n;
3199  return 0;
3200  }
3201  } else { /* decreasing */
3202  if (val > val0) {
3203  *pindex = 0;
3204  return 0;
3205  } else if (val < valn) {
3206  *pindex = n;
3207  return 0;
3208  }
3209  }
3210 
3211  /* Within bounds of existing array; search */
3212  lindex = 0;
3213  rindex = n - 1;
3214  while (1) {
3215  midindex = (lindex + rindex) / 2;
3216  if (midindex == lindex || midindex == rindex) break;
3217  numaGetFValue(na, midindex, &valmid);
3218  if (increasing) {
3219  if (val > valmid)
3220  lindex = midindex;
3221  else
3222  rindex = midindex;
3223  } else { /* decreasing */
3224  if (val > valmid)
3225  rindex = midindex;
3226  else
3227  lindex = midindex;
3228  }
3229  }
3230  *pindex = rindex;
3231  return 0;
3232 }
3233 
3234 
3235 /*----------------------------------------------------------------------*
3236  * Random permutation *
3237  *----------------------------------------------------------------------*/
3253 NUMA *
3255  l_int32 seed)
3256 {
3257 l_int32 i, index, temp;
3258 l_int32 *array;
3259 NUMA *na;
3260 
3261  PROCNAME("numaPseudorandomSequence");
3262 
3263  if (size <= 0)
3264  return (NUMA *)ERROR_PTR("size <= 0", procName, NULL);
3265 
3266  if ((array = (l_int32 *)LEPT_CALLOC(size, sizeof(l_int32))) == NULL)
3267  return (NUMA *)ERROR_PTR("array not made", procName, NULL);
3268  for (i = 0; i < size; i++)
3269  array[i] = i;
3270  srand(seed);
3271  for (i = size - 1; i > 0; i--) {
3272  index = (l_int32)((i + 1) * ((l_float64)rand() / (l_float64)RAND_MAX));
3273  index = L_MIN(index, i);
3274  temp = array[i];
3275  array[i] = array[index];
3276  array[index] = temp;
3277  }
3278 
3279  na = numaCreateFromIArray(array, size);
3280  LEPT_FREE(array);
3281  return na;
3282 }
3283 
3284 
3292 NUMA *
3294  l_int32 seed)
3295 {
3296 l_int32 i, index, size;
3297 l_float32 val;
3298 NUMA *naindex, *nad;
3299 
3300  PROCNAME("numaRandomPermutation");
3301 
3302  if (!nas)
3303  return (NUMA *)ERROR_PTR("nas not defined", procName, NULL);
3304  if ((size = numaGetCount(nas)) == 0) {
3305  L_WARNING("nas is empty\n", procName);
3306  return numaCopy(nas);
3307  }
3308 
3309  naindex = numaPseudorandomSequence(size, seed);
3310  nad = numaCreate(size);
3311  for (i = 0; i < size; i++) {
3312  numaGetIValue(naindex, i, &index);
3313  numaGetFValue(nas, index, &val);
3314  numaAddNumber(nad, val);
3315  }
3316  numaDestroy(&naindex);
3317  return nad;
3318 }
3319 
3320 
3321 /*----------------------------------------------------------------------*
3322  * Functions requiring sorting *
3323  *----------------------------------------------------------------------*/
3351 l_ok
3353  l_float32 fract,
3354  NUMA *nasort,
3355  l_int32 usebins,
3356  l_float32 *pval)
3357 {
3358 l_int32 n, index;
3359 NUMA *nas;
3360 
3361  PROCNAME("numaGetRankValue");
3362 
3363  if (!pval)
3364  return ERROR_INT("&val not defined", procName, 1);
3365  *pval = 0.0; /* init */
3366  if (!na)
3367  return ERROR_INT("na not defined", procName, 1);
3368  if ((n = numaGetCount(na)) == 0)
3369  return ERROR_INT("na empty", procName, 1);
3370  if (fract < 0.0 || fract > 1.0)
3371  return ERROR_INT("fract not in [0.0 ... 1.0]", procName, 1);
3372 
3373  if (nasort) {
3374  nas = nasort;
3375  } else {
3376  if (usebins == 0)
3377  nas = numaSort(NULL, na, L_SORT_INCREASING);
3378  else
3379  nas = numaBinSort(na, L_SORT_INCREASING);
3380  if (!nas)
3381  return ERROR_INT("nas not made", procName, 1);
3382  }
3383  index = (l_int32)(fract * (l_float32)(n - 1) + 0.5);
3384  numaGetFValue(nas, index, pval);
3385 
3386  if (!nasort) numaDestroy(&nas);
3387  return 0;
3388 }
3389 
3390 
3404 l_ok
3406  l_float32 *pval)
3407 {
3408  PROCNAME("numaGetMedian");
3409 
3410  if (!pval)
3411  return ERROR_INT("&val not defined", procName, 1);
3412  *pval = 0.0; /* init */
3413  if (!na || numaGetCount(na) == 0)
3414  return ERROR_INT("na not defined or empty", procName, 1);
3415 
3416  return numaGetRankValue(na, 0.5, NULL, 0, pval);
3417 }
3418 
3419 
3435 l_ok
3437  l_int32 *pval)
3438 {
3439 l_int32 ret;
3440 l_float32 fval;
3441 
3442  PROCNAME("numaGetBinnedMedian");
3443 
3444  if (!pval)
3445  return ERROR_INT("&val not defined", procName, 1);
3446  *pval = 0; /* init */
3447  if (!na || numaGetCount(na) == 0)
3448  return ERROR_INT("na not defined or empty", procName, 1);
3449 
3450  ret = numaGetRankValue(na, 0.5, NULL, 1, &fval);
3451  *pval = lept_roundftoi(fval);
3452  return ret;
3453 }
3454 
3455 
3464 l_ok
3466  l_float32 med,
3467  l_float32 *pdev)
3468 {
3469 l_int32 i, n;
3470 l_float32 val, dev;
3471 
3472  PROCNAME("numaGetMeanDevFromMedian");
3473 
3474  if (!pdev)
3475  return ERROR_INT("&dev not defined", procName, 1);
3476  *pdev = 0.0; /* init */
3477  if (!na)
3478  return ERROR_INT("na not defined", procName, 1);
3479  if ((n = numaGetCount(na)) == 0)
3480  return ERROR_INT("na is empty", procName, 1);
3481 
3482  dev = 0.0;
3483  for (i = 0; i < n; i++) {
3484  numaGetFValue(na, i, &val);
3485  dev += L_ABS(val - med);
3486  }
3487  *pdev = dev / (l_float32)n;
3488  return 0;
3489 }
3490 
3491 
3510 l_ok
3512  l_float32 *pmed,
3513  l_float32 *pdev)
3514 {
3515 l_int32 n, i;
3516 l_float32 val, med;
3517 NUMA *nadev;
3518 
3519  PROCNAME("numaGetMedianDevFromMedian");
3520 
3521  if (pmed) *pmed = 0.0;
3522  if (!pdev)
3523  return ERROR_INT("&dev not defined", procName, 1);
3524  *pdev = 0.0;
3525  if (!na || numaGetCount(na) == 0)
3526  return ERROR_INT("na not defined or empty", procName, 1);
3527 
3528  numaGetMedian(na, &med);
3529  if (pmed) *pmed = med;
3530  n = numaGetCount(na);
3531  nadev = numaCreate(n);
3532  for (i = 0; i < n; i++) {
3533  numaGetFValue(na, i, &val);
3534  numaAddNumber(nadev, L_ABS(val - med));
3535  }
3536  numaGetMedian(nadev, pdev);
3537 
3538  numaDestroy(&nadev);
3539  return 0;
3540 }
3541 
3542 
3559 l_ok
3561  l_float32 *pval,
3562  l_int32 *pcount)
3563 {
3564 l_int32 i, n, maxcount, prevcount;
3565 l_float32 val, maxval, prevval;
3566 l_float32 *array;
3567 NUMA *nasort;
3568 
3569  PROCNAME("numaGetMode");
3570 
3571  if (pcount) *pcount = 0;
3572  if (!pval)
3573  return ERROR_INT("&val not defined", procName, 1);
3574  *pval = 0.0;
3575  if (!na)
3576  return ERROR_INT("na not defined", procName, 1);
3577  if ((n = numaGetCount(na)) == 0)
3578  return ERROR_INT("na is empty", procName, 1);
3579 
3580  if ((nasort = numaSort(NULL, na, L_SORT_DECREASING)) == NULL)
3581  return ERROR_INT("nas not made", procName, 1);
3582  array = numaGetFArray(nasort, L_NOCOPY);
3583 
3584  /* Initialize with array[0] */
3585  prevval = array[0];
3586  prevcount = 1;
3587  maxval = prevval;
3588  maxcount = prevcount;
3589 
3590  /* Scan the sorted array, aggregating duplicates */
3591  for (i = 1; i < n; i++) {
3592  val = array[i];
3593  if (val == prevval) {
3594  prevcount++;
3595  } else { /* new value */
3596  if (prevcount > maxcount) { /* new max */
3597  maxcount = prevcount;
3598  maxval = prevval;
3599  }
3600  prevval = val;
3601  prevcount = 1;
3602  }
3603  }
3604 
3605  /* Was the mode the last run of elements? */
3606  if (prevcount > maxcount) {
3607  maxcount = prevcount;
3608  maxval = prevval;
3609  }
3610 
3611  *pval = maxval;
3612  if (pcount)
3613  *pcount = maxcount;
3614 
3615  numaDestroy(&nasort);
3616  return 0;
3617 }
3618 
3619 
3620 /*----------------------------------------------------------------------*
3621  * Rearrangements *
3622  *----------------------------------------------------------------------*/
3639 l_ok
3641  NUMA *nas,
3642  l_int32 istart,
3643  l_int32 iend)
3644 {
3645 l_int32 n, i;
3646 l_float32 val;
3647 
3648  PROCNAME("numaJoin");
3649 
3650  if (!nad)
3651  return ERROR_INT("nad not defined", procName, 1);
3652  if (!nas)
3653  return 0;
3654 
3655  if (istart < 0)
3656  istart = 0;
3657  n = numaGetCount(nas);
3658  if (iend < 0 || iend >= n)
3659  iend = n - 1;
3660  if (istart > iend)
3661  return ERROR_INT("istart > iend; nothing to add", procName, 1);
3662 
3663  for (i = istart; i <= iend; i++) {
3664  numaGetFValue(nas, i, &val);
3665  numaAddNumber(nad, val);
3666  }
3667 
3668  return 0;
3669 }
3670 
3671 
3688 l_ok
3690  NUMAA *naas,
3691  l_int32 istart,
3692  l_int32 iend)
3693 {
3694 l_int32 n, i;
3695 NUMA *na;
3696 
3697  PROCNAME("numaaJoin");
3698 
3699  if (!naad)
3700  return ERROR_INT("naad not defined", procName, 1);
3701  if (!naas)
3702  return 0;
3703 
3704  if (istart < 0)
3705  istart = 0;
3706  n = numaaGetCount(naas);
3707  if (iend < 0 || iend >= n)
3708  iend = n - 1;
3709  if (istart > iend)
3710  return ERROR_INT("istart > iend; nothing to add", procName, 1);
3711 
3712  for (i = istart; i <= iend; i++) {
3713  na = numaaGetNuma(naas, i, L_CLONE);
3714  numaaAddNuma(naad, na, L_INSERT);
3715  }
3716 
3717  return 0;
3718 }
3719 
3720 
3736 NUMA *
3738 {
3739 l_int32 i, nalloc;
3740 NUMA *na, *nad;
3741 NUMA **array;
3742 
3743  PROCNAME("numaaFlattenToNuma");
3744 
3745  if (!naa)
3746  return (NUMA *)ERROR_PTR("naa not defined", procName, NULL);
3747 
3748  nalloc = naa->nalloc;
3749  array = numaaGetPtrArray(naa);
3750  nad = numaCreate(0);
3751  for (i = 0; i < nalloc; i++) {
3752  na = array[i];
3753  if (!na) continue;
3754  numaJoin(nad, na, 0, -1);
3755  }
3756 
3757  return nad;
3758 }
3759 
l_ok numaGetSum(NUMA *na, l_float32 *psum)
numaGetSum()
Definition: numafunc1.c:538
NUMA * numaInvertMap(NUMA *nas)
numaInvertMap()
Definition: numafunc1.c:3068
NUMA * numaBinSort(NUMA *nas, l_int32 sortorder)
numaBinSort()
Definition: numafunc1.c:2718
void * ptraGetPtrToItem(L_PTRA *pa, l_int32 index)
ptraGetPtrToItem()
Definition: ptra.c:767
l_ok numaGetFValue(NUMA *na, l_int32 index, l_float32 *pval)
numaGetFValue()
Definition: numabasic.c:719
l_ok numaFitMax(NUMA *na, l_float32 *pmaxval, NUMA *naloc, l_float32 *pmaxloc)
numaFitMax()
Definition: numafunc1.c:2162
l_ok numaAddToNumber(NUMA *na, l_int32 index, l_float32 val)
numaAddToNumber()
Definition: numafunc1.c:419
l_ok numaHasOnlyIntegers(NUMA *na, l_int32 *pallints)
numaHasOnlyIntegers()
Definition: numafunc1.c:656
NUMA * numaLowPassIntervals(NUMA *nas, l_float32 thresh, l_float32 maxn)
numaLowPassIntervals()
Definition: numafunc1.c:1410
NUMA * numaGetSortIndex(NUMA *na, l_int32 sortorder)
numaGetSortIndex()
Definition: numafunc1.c:2751
NUMA * numaSortAutoSelect(NUMA *nas, l_int32 sortorder)
numaSortAutoSelect()
Definition: numafunc1.c:2521
l_ok numaInterpolateEqxVal(l_float32 startx, l_float32 deltax, NUMA *nay, l_int32 type, l_float32 xval, l_float32 *pyval)
numaInterpolateEqxVal()
Definition: numafunc1.c:1703
l_int32 lept_roundftoi(l_float32 fval)
lept_roundftoi()
Definition: utils1.c:700
NUMA * numaReverse(NUMA *nad, NUMA *nas)
numaReverse()
Definition: numafunc1.c:1355
NUMA ** numaaGetPtrArray(NUMAA *naa)
numaaGetPtrArray()
Definition: numabasic.c:1719
NUMA * numaSortByIndex(NUMA *nas, NUMA *naindex)
numaSortByIndex()
Definition: numafunc1.c:2916
l_ok numaInsertNumber(NUMA *na, l_int32 index, l_float32 val)
numaInsertNumber()
Definition: numabasic.c:553
Definition: pix.h:713
l_ok numaAddNumber(NUMA *na, l_float32 val)
numaAddNumber()
Definition: numabasic.c:478
NUMA * numaRandomPermutation(NUMA *nas, l_int32 seed)
numaRandomPermutation()
Definition: numafunc1.c:3293
l_ok numaSortPair(NUMA *nax, NUMA *nay, l_int32 sortorder, NUMA **pnasx, NUMA **pnasy)
numaSortPair()
Definition: numafunc1.c:3017
NUMA * numaaFlattenToNuma(NUMAA *naa)
numaaFlattenToNuma()
Definition: numafunc1.c:3737
l_ok numaGetNonzeroRange(NUMA *na, l_float32 eps, l_int32 *pfirst, l_int32 *plast)
numaGetNonzeroRange()
Definition: numafunc1.c:1078
NUMA * numaGetBinSortIndex(NUMA *nas, l_int32 sortorder)
numaGetBinSortIndex()
Definition: numafunc1.c:2833
Definition: pix.h:712
l_int32 numaIsSorted(NUMA *nas, l_int32 sortorder, l_int32 *psorted)
numaIsSorted()
Definition: numafunc1.c:2965
l_int32 nalloc
Definition: array.h:84
l_ok numaGetMedian(NUMA *na, l_float32 *pval)
numaGetMedian()
Definition: numafunc1.c:3405
l_ok numaGetCountRelativeToZero(NUMA *na, l_int32 type, l_int32 *pcount)
numaGetCountRelativeToZero()
Definition: numafunc1.c:1131
NUMA * numaMakeConstant(l_float32 val, l_int32 size)
numaMakeConstant()
Definition: numafunc1.c:851
void lept_stderr(const char *fmt,...)
lept_stderr()
Definition: utils1.c:306
Definition: pix.h:710
l_ok numaIntegrateInterval(NUMA *nax, NUMA *nay, l_float32 x0, l_float32 x1, l_int32 npts, l_float32 *psum)
numaIntegrateInterval()
Definition: numafunc1.c:2351
l_float32 startx
Definition: array.h:75
NUMA * numaSort(NUMA *naout, NUMA *nain, l_int32 sortorder)
numaSort()
Definition: numafunc1.c:2650
l_int32 numaChooseSortType(NUMA *nas)
numaChooseSortType()
Definition: numafunc1.c:2602
l_ok numaCountNonzeroRuns(NUMA *na, l_int32 *pcount)
numaCountNonzeroRuns()
Definition: numafunc1.c:1037
NUMA * numaCreate(l_int32 n)
numaCreate()
Definition: numabasic.c:194
l_float32 delx
Definition: array.h:76
NUMA * numaThresholdEdges(NUMA *nas, l_float32 thresh1, l_float32 thresh2, l_float32 maxn)
numaThresholdEdges()
Definition: numafunc1.c:1487
void * ptraRemoveLast(L_PTRA *pa)
ptraRemoveLast()
Definition: ptra.c:491
l_int32 numaGetSpanValues(NUMA *na, l_int32 span, l_int32 *pstart, l_int32 *pend)
numaGetSpanValues()
Definition: numafunc1.c:1608
NUMA * numaClipToInterval(NUMA *nas, l_int32 first, l_int32 last)
numaClipToInterval()
Definition: numafunc1.c:1182
NUMA * numaGetPartialSums(NUMA *na)
numaGetPartialSums()
Definition: numafunc1.c:579
NUMA * numaSubsample(NUMA *nas, l_int32 subfactor)
numaSubsample()
Definition: numafunc1.c:751
NUMA * numaRemoveBorder(NUMA *nas, l_int32 left, l_int32 right)
numaRemoveBorder()
Definition: numafunc1.c:996
l_ok numaSetValue(NUMA *na, l_int32 index, l_float32 val)
numaSetValue()
Definition: numabasic.c:786
l_ok ptraInsert(L_PTRA *pa, l_int32 index, void *item, l_int32 shiftflag)
ptraInsert()
Definition: ptra.c:344
NUMA * numaUniformSampling(NUMA *nas, l_int32 nsamp)
numaUniformSampling()
Definition: numafunc1.c:1289
NUMA * numaArithOp(NUMA *nad, NUMA *na1, NUMA *na2, l_int32 op)
numaArithOp()
Definition: numafunc1.c:173
NUMA * numaPseudorandomSequence(l_int32 size, l_int32 seed)
numaPseudorandomSequence()
Definition: numafunc1.c:3254
NUMA * numaaGetNuma(NUMAA *naa, l_int32 index, l_int32 accessflag)
numaaGetNuma()
Definition: numabasic.c:1740
NUMA * numaMakeDelta(NUMA *nas)
numaMakeDelta()
Definition: numafunc1.c:786
l_ok numaGetIValue(NUMA *na, l_int32 index, l_int32 *pival)
numaGetIValue()
Definition: numabasic.c:754
Definition: array.h:70
l_ok numaDifferentiateInterval(NUMA *nax, NUMA *nay, l_float32 x0, l_float32 x1, l_int32 npts, NUMA **pnadx, NUMA **pnady)
numaDifferentiateInterval()
Definition: numafunc1.c:2268
L_PTRA * ptraCreate(l_int32 n)
ptraCreate()
Definition: ptra.c:144
l_ok numaInterpolateArbxVal(NUMA *nax, NUMA *nay, l_int32 type, l_float32 xval, l_float32 *pyval)
numaInterpolateArbxVal()
Definition: numafunc1.c:1795
l_int32 numaGetCount(NUMA *na)
numaGetCount()
Definition: numabasic.c:658
l_ok numaGetMeanDevFromMedian(NUMA *na, l_float32 med, l_float32 *pdev)
numaGetMeanDevFromMedian()
Definition: numafunc1.c:3465
l_ok numaGetMode(NUMA *na, l_float32 *pval, l_int32 *pcount)
numaGetMode()
Definition: numafunc1.c:3560
l_ok numaAddSorted(NUMA *na, l_float32 val)
numaAddSorted()
Definition: numafunc1.c:3124
l_ok ptraGetMaxIndex(L_PTRA *pa, l_int32 *pmaxindex)
ptraGetMaxIndex()
Definition: ptra.c:707
NUMA * numaMakeAbsval(NUMA *nad, NUMA *nas)
numaMakeAbsval()
Definition: numafunc1.c:867
l_float32 * array
Definition: array.h:77
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()
Definition: numafunc1.c:1912
Definition: ptra.h:53
l_ok numaFindSortedLoc(NUMA *na, l_float32 val, l_int32 *pindex)
numaFindSortedLoc()
Definition: numafunc1.c:3164
NUMA * numaMakeThresholdIndicator(NUMA *nas, l_float32 thresh, l_int32 type)
numaMakeThresholdIndicator()
Definition: numafunc1.c:1231
l_ok numaSetParameters(NUMA *na, l_float32 startx, l_float32 delx)
numaSetParameters()
Definition: numabasic.c:993
l_int32 numaGetEdgeValues(NUMA *na, l_int32 edge, l_int32 *pstart, l_int32 *pend, l_int32 *psign)
numaGetEdgeValues()
Definition: numafunc1.c:1645
l_ok numaInterpolateArbxInterval(NUMA *nax, NUMA *nay, l_int32 type, l_float32 x0, l_float32 x1, l_int32 npts, NUMA **pnadx, NUMA **pnady)
numaInterpolateArbxInterval()
Definition: numafunc1.c:2001
NUMA * numaAddSpecifiedBorder(NUMA *nas, l_int32 left, l_int32 right, l_int32 type)
numaAddSpecifiedBorder()
Definition: numafunc1.c:945
l_ok numaGetParameters(NUMA *na, l_float32 *pstartx, l_float32 *pdelx)
numaGetParameters()
Definition: numabasic.c:963
Definition: array.h:82
NUMA * numaCreateFromIArray(l_int32 *iarray, l_int32 size)
numaCreateFromIArray()
Definition: numabasic.c:234
l_ok numaGetMean(NUMA *na, l_float32 *pave)
numaGetMean()
Definition: numafunc1.c:691
NUMA * numaMakeSequence(l_float32 startval, l_float32 increment, l_int32 size)
numaMakeSequence()
Definition: numafunc1.c:821
Definition: pix.h:711
l_ok numaSortGeneral(NUMA *na, NUMA **pnasort, NUMA **pnaindex, NUMA **pnainvert, l_int32 sortorder, l_int32 sorttype)
numaSortGeneral()
Definition: numafunc1.c:2456
NUMA * numaCopy(NUMA *na)
numaCopy()
Definition: numabasic.c:399
NUMA * numaAddBorder(NUMA *nas, l_int32 left, l_int32 right, l_float32 val)
numaAddBorder()
Definition: numafunc1.c:902
void numaDestroy(NUMA **pna)
numaDestroy()
Definition: numabasic.c:366
l_ok numaGetMin(NUMA *na, l_float32 *pminval, l_int32 *piminloc)
numaGetMin()
Definition: numafunc1.c:453
void ptraDestroy(L_PTRA **ppa, l_int32 freeflag, l_int32 warnflag)
ptraDestroy()
Definition: ptra.c:194
l_ok numaaJoin(NUMAA *naad, NUMAA *naas, l_int32 istart, l_int32 iend)
numaaJoin()
Definition: numafunc1.c:3689
l_ok numaJoin(NUMA *nad, NUMA *nas, l_int32 istart, l_int32 iend)
numaJoin()
Definition: numafunc1.c:3640
l_float32 * numaGetFArray(NUMA *na, l_int32 copyflag)
numaGetFArray()
Definition: numabasic.c:892
l_int32 numaaGetCount(NUMAA *naa)
numaaGetCount()
Definition: numabasic.c:1631
void * ptraRemove(L_PTRA *pa, l_int32 index, l_int32 flag)
ptraRemove()
Definition: ptra.c:442
NUMA * numaLogicalOp(NUMA *nad, NUMA *na1, NUMA *na2, l_int32 op)
numaLogicalOp()
Definition: numafunc1.c:253
NUMA * numaSortIndexAutoSelect(NUMA *nas, l_int32 sortorder)
numaSortIndexAutoSelect()
Definition: numafunc1.c:2562
l_ok numaGetMedianDevFromMedian(NUMA *na, l_float32 *pmed, l_float32 *pdev)
numaGetMedianDevFromMedian()
Definition: numafunc1.c:3511
l_ok numaReplaceNumber(NUMA *na, l_int32 index, l_float32 val)
numaReplaceNumber()
Definition: numabasic.c:627
NUMA * numaClone(NUMA *na)
numaClone()
Definition: numabasic.c:428
l_int32 numaSimilar(NUMA *na1, NUMA *na2, l_float32 maxdiff, l_int32 *psimilar)
numaSimilar()
Definition: numafunc1.c:370
l_ok numaGetMax(NUMA *na, l_float32 *pmaxval, l_int32 *pimaxloc)
numaGetMax()
Definition: numafunc1.c:496
l_ok numaGetBinnedMedian(NUMA *na, l_int32 *pval)
numaGetBinnedMedian()
Definition: numafunc1.c:3436
NUMA * numaInvert(NUMA *nad, NUMA *nas)
numaInvert()
Definition: numafunc1.c:326
l_ok numaGetMeanAbsval(NUMA *na, l_float32 *paveabs)
numaGetMeanAbsval()
Definition: numafunc1.c:720
l_ok numaGetSumOnInterval(NUMA *na, l_int32 first, l_int32 last, l_float32 *psum)
numaGetSumOnInterval()
Definition: numafunc1.c:613
l_ok numaGetRankValue(NUMA *na, l_float32 fract, NUMA *nasort, l_int32 usebins, l_float32 *pval)
numaGetRankValue()
Definition: numafunc1.c:3352
l_ok numaaAddNuma(NUMAA *naa, NUMA *na, l_int32 copyflag)
numaaAddNuma()
Definition: numabasic.c:1546