quant.c 21 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788
  1. /*!
  2. * \file lib/raster/quant.c
  3. *
  4. * \brief Raster Library - Quantization rules.
  5. *
  6. * The quantization table is stored as a linear array. Rules are added
  7. * starting from index 0. Redundant rules are not eliminated. Rules
  8. * are tested from the highest index downto 0. There are two
  9. * "infinite" rules. Support is provided to reverse the order of the
  10. * rules.
  11. *
  12. * (C) 1999-2009 by the GRASS Development Team
  13. *
  14. * This program is free software under the GNU General Public License
  15. * (>=v2). Read the file COPYING that comes with GRASS for details.
  16. *
  17. * \author USACERL and many others
  18. */
  19. #include <stdlib.h>
  20. #include <grass/gis.h>
  21. #include <grass/raster.h>
  22. static int double_comp(const void *, const void *);
  23. #define USE_LOOKUP 1
  24. #define MAX_LOOKUP_TABLE_SIZE 2048
  25. #define NO_DATA (Rast_set_c_null_value (&tmp, 1), (CELL) tmp)
  26. #define NO_LEFT_INFINITE_RULE (! q->infiniteLeftSet)
  27. #define NO_RIGHT_INFINITE_RULE (! q->infiniteRightSet)
  28. #define NO_FINITE_RULE (q->nofRules <= 0)
  29. #define NO_EXPLICIT_RULE (NO_FINITE_RULE && \
  30. NO_LEFT_INFINITE_RULE && NO_RIGHT_INFINITE_RULE)
  31. /*!
  32. \brief Resets the number of defined rules and number of infinite rules to 0
  33. \param q pointer to Quant structure to be reset
  34. */
  35. void Rast_quant_clear(struct Quant *q)
  36. {
  37. q->nofRules = 0;
  38. q->infiniteRightSet = q->infiniteLeftSet = 0;
  39. }
  40. /*!
  41. \brief Resets and frees allocated memory
  42. Resets the number of defined rules to 0 and free's space allocated
  43. for rules. Calls Rast_quant_clear().
  44. \param q pointer to Quant structure to be reset
  45. */
  46. void Rast_quant_free(struct Quant *q)
  47. {
  48. Rast_quant_clear(q);
  49. if (q->maxNofRules > 0)
  50. G_free(q->table);
  51. if (q->fp_lookup.active) {
  52. G_free(q->fp_lookup.vals);
  53. G_free(q->fp_lookup.rules);
  54. q->fp_lookup.nalloc = 0;
  55. q->fp_lookup.active = 0;
  56. }
  57. q->maxNofRules = 0;
  58. }
  59. /*!
  60. * \brief Organized fp_lookup table.
  61. *
  62. * Organizes fp_lookup table for faster (logarithmic) lookup time
  63. * G_quant_organize_fp_lookup() creates a list of min and max for
  64. * each quant rule, sorts this list, and stores the pointer to quant
  65. * rule that should be used inbetween any 2 numbers in this list.
  66. * Also it stores extreme points for 2 infinite rules, if exist.
  67. * After the call to G_quant_organize_fp_lookup()
  68. * instead of linearly searching through list of rules to find
  69. * a rule to apply, quant lookup will perform a binary search
  70. * to find an interval containing floating point value, and then use
  71. * the rule associated with this interval.
  72. * when the value doesn't fall within any interval, check for the
  73. * infinite rules.
  74. *
  75. * \param q pointer to Quant structure which holds quant rules info
  76. *
  77. * \return 1 on success
  78. */
  79. int Rast__quant_organize_fp_lookup(struct Quant *q)
  80. {
  81. int i;
  82. DCELL val;
  83. CELL tmp;
  84. struct Quant_table *p;
  85. if (q->nofRules * 2 > MAX_LOOKUP_TABLE_SIZE)
  86. return -1;
  87. if (q->nofRules == 0)
  88. return -1;
  89. q->fp_lookup.vals = (DCELL *)
  90. G_calloc(q->nofRules * 2, sizeof(DCELL));
  91. /* 2 endpoints for each rule */
  92. q->fp_lookup.rules = (struct Quant_table **)
  93. G_calloc(q->nofRules * 2, sizeof(struct Quant_table *));
  94. /* first we organize finite rules into a table */
  95. if (!NO_FINITE_RULE) {
  96. i = 0;
  97. /* get the list of DCELL values from set of all dLows and dHighs
  98. of all rules */
  99. /* NOTE: if dLow==DHigh in a rule, the value appears twice in a list
  100. but if dLow==DHigh of the previous, rule the value appears only once */
  101. for (p = &(q->table[q->nofRules - 1]); p >= q->table; p--) {
  102. /* check if the min is the same as previous maximum */
  103. if (i == 0 || p->dLow != q->fp_lookup.vals[i - 1])
  104. q->fp_lookup.vals[i++] = p->dLow;
  105. q->fp_lookup.vals[i++] = p->dHigh;
  106. }
  107. q->fp_lookup.nalloc = i;
  108. /* now sort the values */
  109. qsort((char *)q->fp_lookup.vals, q->fp_lookup.nalloc,
  110. sizeof(DCELL), double_comp);
  111. /* now find the rule to apply inbetween each 2 values in a list */
  112. for (i = 0; i < q->fp_lookup.nalloc - 1; i++) {
  113. /*debug
  114. fprintf (stderr, "%lf %lf ", q->fp_lookup.vals[i], q->fp_lookup.vals[i+1]);
  115. */
  116. val = (q->fp_lookup.vals[i] + q->fp_lookup.vals[i + 1]) / 2.;
  117. q->fp_lookup.rules[i] =
  118. Rast__quant_get_rule_for_d_raster_val(q, val);
  119. /* debug
  120. if(q->fp_lookup.rules[i])
  121. fprintf (stderr, "%lf %lf %d %d\n", q->fp_lookup.rules[i]->dLow, q->fp_lookup.rules[i]->dHigh, q->fp_lookup.rules[i]->cLow, q->fp_lookup.rules[i]->cHigh);
  122. else fprintf (stderr, "null\n");
  123. */
  124. }
  125. } /* organizing finite rules */
  126. if (!NO_LEFT_INFINITE_RULE) {
  127. q->fp_lookup.inf_dmin = q->infiniteDLeft;
  128. q->fp_lookup.inf_min = q->infiniteCLeft;
  129. }
  130. else {
  131. if (q->fp_lookup.nalloc)
  132. q->fp_lookup.inf_dmin = q->fp_lookup.vals[0];
  133. q->fp_lookup.inf_min = NO_DATA;
  134. }
  135. if (!NO_RIGHT_INFINITE_RULE) {
  136. if (q->fp_lookup.nalloc)
  137. q->fp_lookup.inf_dmax = q->infiniteDRight;
  138. q->fp_lookup.inf_max = q->infiniteCRight;
  139. }
  140. else {
  141. q->fp_lookup.inf_dmax = q->fp_lookup.vals[q->fp_lookup.nalloc - 1];
  142. q->fp_lookup.inf_max = NO_DATA;
  143. }
  144. q->fp_lookup.active = 1;
  145. return 1;
  146. }
  147. /*!
  148. * \brief Initialize the structure
  149. *
  150. * Initializes the <i>q</i> struct.
  151. *
  152. * \param quant pointer to Quant structure to be initialized
  153. */
  154. void Rast_quant_init(struct Quant *quant)
  155. {
  156. quant->fp_lookup.active = 0;
  157. quant->maxNofRules = 0;
  158. quant->truncate_only = 0;
  159. quant->round_only = 0;
  160. Rast_quant_clear(quant);
  161. }
  162. /*!
  163. \brief Returns whether or not quant rules are set to truncate map
  164. \param quant pointer to Quant structure which holds quant rules info
  165. \return 1 if truncate is enable
  166. \return 0 if not truncated
  167. */
  168. int Rast_quant_is_truncate(const struct Quant *quant)
  169. {
  170. return quant->truncate_only;
  171. }
  172. /*!
  173. \brief Returns whether or not quant rules are set to round map
  174. \param quant pointer to Quant structure which holds quant rules info
  175. \return 1 is round
  176. \return 0 not round
  177. */
  178. int Rast_quant_is_round(const struct Quant *quant)
  179. {
  180. return quant->round_only;
  181. }
  182. /*!
  183. * \brief Sets the quant rules to perform simple truncation on floats.
  184. *
  185. * Sets the quant for <i>q</i> rules to perform simple truncation on
  186. * floats.
  187. *
  188. * \param quant pointer to Quant structure which holds quant rules info
  189. */
  190. void Rast_quant_truncate(struct Quant *quant)
  191. {
  192. quant->truncate_only = 1;
  193. }
  194. /*!
  195. * \brief Sets the quant rules to perform simple rounding on floats.
  196. *
  197. * Sets the quant for <i>q</i> rules to perform simple rounding on
  198. * floats.
  199. *
  200. * \param quant pointer to Quant structure which holds quant rules info
  201. */
  202. void Rast_quant_round(struct Quant *quant)
  203. {
  204. quant->round_only = 1;
  205. }
  206. static void quant_set_limits(struct Quant *q,
  207. DCELL dLow, DCELL dHigh, CELL cLow, CELL cHigh)
  208. {
  209. q->dMin = dLow;
  210. q->dMax = dHigh;
  211. q->cMin = cLow;
  212. q->cMax = cHigh;
  213. }
  214. static void quant_update_limits(struct Quant *q,
  215. DCELL dLow, DCELL dHigh,
  216. CELL cLow, DCELL cHigh)
  217. {
  218. if (NO_EXPLICIT_RULE) {
  219. quant_set_limits(q, dLow, dHigh, cLow, cHigh);
  220. return;
  221. }
  222. q->dMin = MIN(q->dMin, MIN(dLow, dHigh));
  223. q->dMax = MAX(q->dMax, MAX(dLow, dHigh));
  224. q->cMin = MIN(q->cMin, MIN(cLow, cHigh));
  225. q->cMax = MAX(q->cMax, MAX(cLow, cHigh));
  226. }
  227. /*!
  228. * \brief Returns the minimum and maximum cell and dcell values of all
  229. * the ranges defined.
  230. *
  231. * Extracts the minimum and maximum floating-point and integer values
  232. * from all the rules (except the "infinite" rules) in <i>q</i> into
  233. * <i>dmin</i>, <i>dmax</i>, <i>cmin</i>, and <i>cmax</i>.
  234. *
  235. * \param quant pointer to Quant structure which holds quant rules info
  236. * \param[out] dmin minimum fp value
  237. * \param[out] dmax maximum fp value
  238. * \param[out] cmin minimum value
  239. * \param[out] cmax maximum value
  240. *
  241. * \return -1 if q->truncate or q->round are true or after
  242. * Rast_quant_init (), or any call to Rast_quant_clear () or Rast_quant_free()
  243. * no explicit rules have been added. In this case the returned
  244. * minimum and maximum CELL and DCELL values are null.
  245. * \return 1 if there are any explicit rules
  246. * \return 0 if there are no explicit rules (this includes cases when
  247. * q is set to truncate or round map), and sets <i>dmin</i>,
  248. * <i>dmax</i>, <i>cmin</i>, and <i>cmax</i> to NULL.
  249. */
  250. int Rast_quant_get_limits(const struct Quant *q,
  251. DCELL * dMin, DCELL * dMax, CELL * cMin,
  252. CELL * cMax)
  253. {
  254. if (NO_EXPLICIT_RULE) {
  255. Rast_set_c_null_value(cMin, 1);
  256. Rast_set_c_null_value(cMax, 1);
  257. Rast_set_d_null_value(dMin, 1);
  258. Rast_set_d_null_value(dMax, 1);
  259. return -1;
  260. }
  261. *dMin = q->dMin;
  262. *dMax = q->dMax;
  263. *cMin = q->cMin;
  264. *cMax = q->cMax;
  265. return 1;
  266. }
  267. /*!
  268. \brief Returns the number of quantization rules defined.
  269. This number does not include the 2 infinite intervals.
  270. \param q pointer to Quant structure which holds quant rules info
  271. \return number of quantization rules
  272. */
  273. int Rast_quant_nof_rules(const struct Quant *q)
  274. {
  275. return q->nofRules;
  276. }
  277. /*!
  278. \brief Returns the i'th quantization rule.
  279. For 0 <= i < Rast_quant_nof_rules(). A larger value for i means that
  280. the rule has been added later.
  281. \param q pointer to Quant structure which holds quant rules info
  282. \param i index
  283. \param[out] dLow minimum fp value
  284. \param[out] dHigh maximum fp value
  285. \param[out] cLow minimum value
  286. \param[out] cHigh maximum value
  287. */
  288. void Rast_quant_get_ith_rule(const struct Quant *q,
  289. int i,
  290. DCELL * dLow, DCELL * dHigh,
  291. CELL * cLow, CELL * cHigh)
  292. {
  293. *dLow = q->table[i].dLow;
  294. *dHigh = q->table[i].dHigh;
  295. *cLow = q->table[i].cLow;
  296. *cHigh = q->table[i].cHigh;
  297. }
  298. static void quant_table_increase(struct Quant *q)
  299. {
  300. if (q->nofRules < q->maxNofRules)
  301. return;
  302. if (q->maxNofRules == 0) {
  303. q->maxNofRules = 50;
  304. q->table = (struct Quant_table *)
  305. G_malloc(q->maxNofRules * sizeof(struct Quant_table));
  306. }
  307. else {
  308. q->maxNofRules += 50;
  309. q->table = (struct Quant_table *)
  310. G_realloc((char *)q->table,
  311. q->maxNofRules * sizeof(struct Quant_table));
  312. }
  313. }
  314. /*!
  315. \brief Defines a rule for values "dLeft" and smaller.
  316. Values in this range are mapped to "c" if none of the "finite"
  317. quantization rules applies.
  318. \param q pointer to Quant structure which holds quant rules info
  319. \param dLeft fp value
  320. \param c value
  321. */
  322. void Rast_quant_set_neg_infinite_rule(struct Quant *q, DCELL dLeft, CELL c)
  323. {
  324. q->infiniteDLeft = dLeft;
  325. q->infiniteCLeft = c;
  326. quant_update_limits(q, dLeft, dLeft, c, c);
  327. /* update lookup table */
  328. if (q->fp_lookup.active) {
  329. q->fp_lookup.inf_dmin = q->infiniteDLeft;
  330. q->fp_lookup.inf_min = q->infiniteCLeft;
  331. }
  332. q->infiniteLeftSet = 1;
  333. }
  334. /*!
  335. \brief Returns in "dLeft" and "c" the rule values.
  336. For the negative infinite interval (see Rast_quant_set_neg_infinite_rule()).
  337. \param q pointer to Quant structure which holds quant rules info
  338. \param[out] dLeft fp value
  339. \param[out] c value
  340. \return 0 if this rule is not defined
  341. \return 1 otherwise
  342. */
  343. int Rast_quant_get_neg_infinite_rule(const struct Quant *q,
  344. DCELL * dLeft, CELL * c)
  345. {
  346. if (q->infiniteLeftSet == 0)
  347. return 0;
  348. *dLeft = q->infiniteDLeft;
  349. *c = q->infiniteCLeft;
  350. return 1;
  351. }
  352. /*!
  353. \brief Defines a rule for values "dRight" and larger.
  354. Values in this range are mapped to "c" if none of the "finite"
  355. quantization rules or the negative infinite rule applies.
  356. \param q pointer to Quant structure which holds quant rules info
  357. \param dRight fp value
  358. \param c value
  359. */
  360. void Rast_quant_set_pos_infinite_rule(struct Quant *q, DCELL dRight, CELL c)
  361. {
  362. q->infiniteDRight = dRight;
  363. q->infiniteCRight = c;
  364. quant_update_limits(q, dRight, dRight, c, c);
  365. /* update lookup table */
  366. if (q->fp_lookup.active) {
  367. q->fp_lookup.inf_dmax = q->infiniteDRight;
  368. q->fp_lookup.inf_max = q->infiniteCRight;
  369. }
  370. q->infiniteRightSet = 1;
  371. }
  372. /*!
  373. \brief Returns in "dRight" and "c" the rule values.
  374. For the positive infinite interval (see Rast_quant_set_pos_infinite_rule()).
  375. \param q pointer to Quant structure which holds quant rules info
  376. \param[out] dRight fp value
  377. \param[out] c value
  378. \return 0 if this rule is not defined
  379. \return 1 otherwise
  380. */
  381. int Rast_quant_get_pos_infinite_rule(const struct Quant *q,
  382. DCELL * dRight, CELL * c)
  383. {
  384. if (q->infiniteRightSet == 0)
  385. return 0;
  386. *dRight = q->infiniteDRight;
  387. *c = q->infiniteCRight;
  388. return 1;
  389. }
  390. /*!
  391. \brief Adds a new rule to the set of quantization rules.
  392. If dLow < dHigh the rule will be stored with the low and high values
  393. interchanged.
  394. Note: currently no cleanup of rules is performed, i.e. redundant
  395. rules are not removed. This can't be changed because Categories
  396. structure HEAVILY depends of quant rules stored in exactly the same
  397. order they are entered. So if the cleanup or rearrangement is done in
  398. the future make a flag for add_rule whether or not to do it, then
  399. quant will not set this flag.
  400. \param q pointer to Quant structure which holds quant rules info
  401. \param dLow minimum fp value
  402. \param dHigh maximum fp value
  403. \param cLow minimum value
  404. \param cHigh maximum value
  405. */
  406. void Rast_quant_add_rule(struct Quant *q,
  407. DCELL dLow, DCELL dHigh, CELL cLow, CELL cHigh)
  408. {
  409. int i;
  410. struct Quant_table *p;
  411. quant_table_increase(q);
  412. i = q->nofRules;
  413. p = &(q->table[i]);
  414. if (dHigh >= dLow) {
  415. p->dLow = dLow;
  416. p->dHigh = dHigh;
  417. p->cLow = cLow;
  418. p->cHigh = cHigh;
  419. }
  420. else {
  421. p->dLow = dHigh;
  422. p->dHigh = dLow;
  423. p->cLow = cHigh;
  424. p->cHigh = cLow;
  425. }
  426. /* destroy lookup table, it has to be rebuilt */
  427. if (q->fp_lookup.active) {
  428. G_free(q->fp_lookup.vals);
  429. G_free(q->fp_lookup.rules);
  430. q->fp_lookup.active = 0;
  431. q->fp_lookup.nalloc = 0;
  432. }
  433. quant_update_limits(q, dLow, dHigh, cLow, cHigh);
  434. q->nofRules++;
  435. }
  436. /*!
  437. \brief Rreverses the order in which the qunatization rules are stored.
  438. See also Rast_quant_get_ith_rule() and Rast_quant_perform_d()).
  439. \param q pointer to Quant rules which holds quant rules info
  440. */
  441. void Rast_quant_reverse_rule_order(struct Quant *q)
  442. {
  443. struct Quant_table tmp;
  444. struct Quant_table *pLeft, *pRight;
  445. pLeft = q->table;
  446. pRight = &(q->table[q->nofRules - 1]);
  447. while (pLeft < pRight) {
  448. tmp.dLow = pLeft->dLow;
  449. tmp.dHigh = pLeft->dHigh;
  450. tmp.cLow = pLeft->cLow;
  451. tmp.cHigh = pLeft->cHigh;
  452. pLeft->dLow = pRight->dLow;
  453. pLeft->dHigh = pRight->dHigh;
  454. pLeft->cLow = pRight->cLow;
  455. pLeft->cHigh = pRight->cHigh;
  456. pRight->dLow = tmp.dLow;
  457. pRight->dHigh = tmp.dHigh;
  458. pRight->cLow = tmp.cLow;
  459. pRight->cHigh = tmp.cHigh;
  460. pLeft++;
  461. pRight--;
  462. }
  463. }
  464. static CELL quant_interpolate(DCELL dLow, DCELL dHigh,
  465. CELL cLow, CELL cHigh, DCELL dValue)
  466. {
  467. if (cLow == cHigh)
  468. return cLow;
  469. if (dLow == dHigh)
  470. return cLow;
  471. return (CELL) ((dValue - dLow) / (dHigh - dLow) * (DCELL) (cHigh - cLow) +
  472. (DCELL) cLow);
  473. }
  474. static int less_or_equal(double x, double y)
  475. {
  476. if (x <= y)
  477. return 1;
  478. else
  479. return 0;
  480. }
  481. static int less(double x, double y)
  482. {
  483. if (x < y)
  484. return 1;
  485. else
  486. return 0;
  487. }
  488. /*!
  489. * \brief
  490. *
  491. *
  492. * Returns a CELL category for the floating-point <i>value</i> based
  493. * on the quantization rules in <i>q</i>. The first rule found that
  494. * applies is used. The rules are searched in the reverse order they
  495. * are added to <i>q</i>. If no rule is found, the <i>value</i>
  496. * is first tested against the negative infinite rule, and finally
  497. * against the positive infinite rule. If none of these rules apply,
  498. * the NULL-value is returned.
  499. *
  500. * <b>Note:</b> See G_quant_organize_fp_lookup() for details on how
  501. * the values are looked up from fp_lookup table when it is
  502. * active. Right now fp_lookup is automatically organized during the
  503. * first call to Rast_quant_get_cell_value().
  504. *
  505. * \param q pointer to Quant structure which holds quant rules info
  506. * \param dcellValue fp cell value
  507. *
  508. * \return cell value (integer)
  509. */
  510. CELL Rast_quant_get_cell_value(struct Quant * q, DCELL dcellVal)
  511. {
  512. CELL tmp;
  513. DCELL dtmp;
  514. int try, min_ind, max_ind;
  515. struct Quant_table *p;
  516. int (*lower) ();
  517. dtmp = dcellVal;
  518. /* I know the functions which call me already check for null values,
  519. but I am a public function, and can be called from outside */
  520. if (Rast_is_d_null_value(&dtmp))
  521. return NO_DATA;
  522. if (q->truncate_only)
  523. return (CELL) dtmp;
  524. if (q->round_only) {
  525. if (dcellVal > 0)
  526. return (CELL) (dcellVal + .5);
  527. return (CELL) (dcellVal - .5);
  528. }
  529. if (NO_EXPLICIT_RULE)
  530. return NO_DATA;
  531. if (NO_EXPLICIT_RULE)
  532. return NO_DATA;
  533. if (USE_LOOKUP &&
  534. (q->fp_lookup.active || Rast__quant_organize_fp_lookup(q) > 0)) {
  535. /* first check if values fall within range */
  536. /* if value is below the range */
  537. if (dcellVal < q->fp_lookup.vals[0]) {
  538. if (dcellVal <= q->fp_lookup.inf_dmin)
  539. return q->fp_lookup.inf_min;
  540. else
  541. return NO_DATA;
  542. }
  543. /* if value is below above range */
  544. if (dcellVal > q->fp_lookup.vals[q->fp_lookup.nalloc - 1]) {
  545. if (dcellVal >= q->fp_lookup.inf_dmax)
  546. return q->fp_lookup.inf_max;
  547. else
  548. return NO_DATA;
  549. }
  550. /* make binary search to find which interval our value belongs to
  551. and apply the rule for this interval */
  552. try = (q->fp_lookup.nalloc - 1) / 2;
  553. min_ind = 0;
  554. max_ind = q->fp_lookup.nalloc - 2;
  555. while (1) {
  556. /* DEBUG
  557. fprintf (stderr, "%d %d %d\n", min_ind, max_ind, try);
  558. */
  559. /* when the ruke for the interval is NULL, we exclude the end points.
  560. when it exists, we include the end-points */
  561. if (q->fp_lookup.rules[try])
  562. lower = less;
  563. else
  564. lower = less_or_equal;
  565. if (lower(q->fp_lookup.vals[try + 1], dcellVal)) { /* recurse to the second half */
  566. min_ind = try + 1;
  567. /* must be still < nalloc-1, since number is within the range */
  568. try = (max_ind + min_ind) / 2;
  569. continue;
  570. }
  571. if (lower(dcellVal, q->fp_lookup.vals[try])) { /* recurse to the second half */
  572. max_ind = try - 1;
  573. /* must be still >= 0, since number is within the range */
  574. try = (max_ind + min_ind) / 2;
  575. continue;
  576. }
  577. /* the value fits into the interval! */
  578. p = q->fp_lookup.rules[try];
  579. if (p)
  580. return quant_interpolate(p->dLow, p->dHigh, p->cLow, p->cHigh,
  581. dcellVal);
  582. /* otherwise when finite rule for this interval doesn't exist */
  583. else { /* first check if maybe infinite rule applies */
  584. if (dcellVal <= q->fp_lookup.inf_dmin)
  585. return q->fp_lookup.inf_min;
  586. if (dcellVal >= q->fp_lookup.inf_dmax)
  587. return q->fp_lookup.inf_max;
  588. else
  589. return NO_DATA;
  590. }
  591. } /* while */
  592. } /* looking up in fp_lookup */
  593. if (!NO_FINITE_RULE) {
  594. p = Rast__quant_get_rule_for_d_raster_val(q, dcellVal);
  595. if (!p)
  596. return NO_DATA;
  597. return quant_interpolate(p->dLow, p->dHigh, p->cLow, p->cHigh,
  598. dcellVal);
  599. }
  600. if ((!NO_LEFT_INFINITE_RULE) && (dcellVal <= q->infiniteDLeft))
  601. return q->infiniteCLeft;
  602. if ((NO_RIGHT_INFINITE_RULE) || (dcellVal < q->infiniteDRight))
  603. return NO_DATA;
  604. return q->infiniteCRight;
  605. }
  606. /*!
  607. \brief Returns in "cell" the quantized CELL values.
  608. Returns in "cell" the quantized CELL values corresponding to the
  609. DCELL values stored in "dcell". the number of elements quantized
  610. is n. quantization is performed by repeated application of
  611. Rast_quant_get_cell_value().
  612. \param q pointer to Quant structure which holds quant rules info
  613. \param dcell pointer to fp cell values array
  614. \param[out] cell pointer cell values array
  615. \param n number of cells
  616. */
  617. void Rast_quant_perform_d(struct Quant *q,
  618. const DCELL * dcell, CELL * cell, int n)
  619. {
  620. int i;
  621. for (i = 0; i < n; i++, dcell++)
  622. if (!Rast_is_d_null_value(dcell))
  623. *cell++ = Rast_quant_get_cell_value(q, *dcell);
  624. else
  625. Rast_set_c_null_value(cell++, 1);
  626. }
  627. /*!
  628. \brief Same as Rast_quant_perform_d(), except the type.
  629. \param q pointer to Quant structure which holds quant rules info
  630. \param fcell pointer to fp cell values array
  631. \param[out] cell pointer cell values array
  632. \param n number of cells
  633. */
  634. void Rast_quant_perform_f(struct Quant *q,
  635. const FCELL * fcell, CELL * cell, int n)
  636. {
  637. int i;
  638. for (i = 0; i < n; i++, fcell++)
  639. if (!Rast_is_f_null_value(fcell))
  640. *cell++ = Rast_quant_get_cell_value(q, (DCELL) * fcell);
  641. else
  642. Rast_set_c_null_value(cell++, 1);
  643. }
  644. static int double_comp(const void *xx, const void *yy)
  645. {
  646. const DCELL *x = xx;
  647. const DCELL *y = yy;
  648. if (Rast_is_d_null_value(x))
  649. return 0;
  650. if (*x < *y)
  651. return -1;
  652. else if (*x == *y)
  653. return 0;
  654. else
  655. return 1;
  656. }
  657. /*!
  658. \brief Returns quant rule which will be applied.
  659. Returns quant rule which will be applied when looking up the integer
  660. quant value for val (used when organizing fp_lookup).
  661. \param q pointer to Quant structure which holds quant rules info
  662. \param val fp cell value
  663. \return pointer to the Quant_table (color rule)
  664. \return NULL otherwise
  665. */
  666. struct Quant_table *Rast__quant_get_rule_for_d_raster_val(const struct Quant
  667. *q, DCELL val)
  668. {
  669. const struct Quant_table *p;
  670. for (p = &(q->table[q->nofRules - 1]); p >= q->table; p--)
  671. if ((val >= p->dLow) && (val <= p->dHigh))
  672. break;
  673. if (p >= q->table)
  674. return (struct Quant_table *)p;
  675. else
  676. return (struct Quant_table *)NULL;
  677. }