OmniSciDB  1dac507f6e
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
ExtensionFunctions.hpp
Go to the documentation of this file.
1 #include "../Shared/funcannotations.h"
2 #ifndef __CUDACC__
3 #include <cstdint>
4 #endif
5 #include <cmath>
6 #include <cstdlib>
7 
8 #define EXTENSION_INLINE extern "C" ALWAYS_INLINE DEVICE
9 #define EXTENSION_NOINLINE extern "C" NEVER_INLINE DEVICE
10 
11 /* Example extension functions:
12  *
13  * EXTENSION_NOINLINE
14  * int32_t diff(const int32_t x, const int32_t y) {
15  * return x - y;
16  *
17  * Arrays map to a pair of pointer and element count. The pointer type must be
18  * consistent with the element type of the array column. Only array of numbers
19  * are supported for the moment. For example, an ARRAY_AT function which works
20  * for arrays of INTEGER and SMALLINT can be implemented as follows:
21  *
22  * EXTENSION_NOINLINE
23  * int64_t array_at(const int32_t* arr, const size_t elem_count, const size_t idx) {
24  * return idx < elem_count ? arr[idx] : -1;
25  * }
26  *
27  * EXTENSION_NOINLINE
28  * int64_t array_at__(const int16_t* arr, const size_t elem_count, const size_t idx) {
29  * return idx < elem_count ? arr[idx] : -1;
30  * }
31  *
32  * Note that the return type cannot be specialized and must be the same across
33  * all specializations. We can remove this constraint in the future if necessary.
34  */
35 
37 double Acos(const double x) {
38  return acos(x);
39 }
40 
42 double Asin(const double x) {
43  return asin(x);
44 }
45 
47 double Atan(const double x) {
48  return atan(x);
49 }
50 
52 double Atan2(const double y, const double x) {
53  return atan2(y, x);
54 }
55 
57 double Ceil(double x) {
58  return ceil(x);
59 }
60 
62 float Ceil__(float x) {
63  return ceil(x);
64 }
65 
67 int16_t Ceil__1(int16_t x) {
68  return x;
69 }
70 
72 int32_t Ceil__2(int32_t x) {
73  return x;
74 }
75 
77 int64_t Ceil__3(int64_t x) {
78  return x;
79 }
80 
82 double Cos(const double x) {
83  return cos(x);
84 }
85 
87 double Cot(const double x) {
88  return 1 / tan(x);
89 }
90 
92 double degrees(double x) {
93  return x * (180.0 / M_PI);
94 }
95 
97 double Exp(double x) {
98  return exp(x);
99 }
100 
102 double Floor(double x) {
103  return floor(x);
104 }
105 
107 float Floor__(float x) {
108  return floor(x);
109 }
110 
112 int16_t Floor__1(int16_t x) {
113  return x;
114 }
115 
117 int32_t Floor__2(int32_t x) {
118  return x;
119 }
120 
122 int64_t Floor__3(int64_t x) {
123  return x;
124 }
125 
127 double ln(const double x) {
128  return log(x);
129 }
130 
132 double ln__(const float x) {
133  return logf(x);
134 }
135 
137 double Log(const double x) {
138  return log(x);
139 }
140 
142 double Log__(const float x) {
143  return logf(x);
144 }
145 
147 double Log10(const double x) {
148  return log10(x);
149 }
150 
152 double Log10__(const float x) {
153  return log10f(x);
154 }
155 
157 double pi() {
158  return M_PI;
159 }
160 
162 double power(const double x, const double y) {
163  return pow(x, y);
164 }
165 
167 double radians(const double x) {
168  return x * (M_PI / 180.0);
169 }
170 
172 double Round(const double x, const int32_t y) {
173  if (y == 0) {
174  return round(x) + 0.0;
175  }
176 
177  double exp = pow(10, y);
178 #if defined(__powerpc__) && !defined(__CUDACC__)
179  int32_t yy = y - 1;
180  exp = 10 * powf((float)10L, yy);
181 #endif
182  return (round(x * exp) / exp) + 0.0;
183 }
184 
186 float Round__(const float x, const int32_t y) {
187  if (y == 0) {
188  return roundf(x) + 0.0f;
189  }
190 
191  float exp = powf((float)10L, y);
192 #if defined(__powerpc__) && !defined(__CUDACC__)
193  int32_t yy = y - 1;
194  exp = 10 * powf((float)10L, yy);
195 #endif
196  return roundf(x * exp) / exp + 0.0f;
197 }
198 
200 int16_t Round__1(const int16_t x, const int32_t y) {
201  if (y >= 0) {
202  return x;
203  }
204 
205  int32_t p = pow((float)10L, std::abs(y));
206  int32_t p_half = p >> 1;
207 
208  int64_t temp = x;
209 #if defined(__powerpc__) && !defined(__CUDACC__)
210  int16_t xx = x;
211  xx += 1;
212  temp = xx;
213  temp -= 1;
214 #endif
215  temp = temp >= 0 ? temp + p_half : temp - p_half;
216  temp = temp / p;
217  return temp * p;
218 }
219 
221 int32_t Round__2(const int32_t x, const int32_t y) {
222  if (y >= 0) {
223  return x;
224  }
225 
226  int32_t p = pow((float)10L, std::abs(y));
227  int32_t p_half = p >> 1;
228 
229  int64_t temp = x;
230 #if defined(__powerpc__) && !defined(__CUDACC__)
231  int32_t xx = x;
232  xx += 1;
233  temp = xx;
234  temp -= 1;
235 #endif
236  temp = temp >= 0 ? temp + p_half : temp - p_half;
237  temp = temp / p;
238  return temp * p;
239 }
240 
242 int64_t Round__3(const int64_t x, const int32_t y) {
243  if (y >= 0) {
244  return x;
245  }
246 
247  int64_t p = pow((double)10L, std::abs(y));
248  int64_t p_half = p >> 1;
249 
250  int64_t temp = x;
251  temp = temp >= 0 ? temp + p_half : temp - p_half;
252  temp = temp / p;
253  return temp * p;
254 }
255 
257 int64_t Round__4(const int64_t x, const int32_t y0, const int32_t scale) {
258  int32_t y = y0 - scale;
259 
260  if (y >= 0) {
261  return x;
262  }
263 
264  int64_t p = pow((double)10L, std::abs(y));
265  int64_t p_half = p >> 1;
266 
267  int64_t temp = x;
268  temp = temp >= 0 ? temp + p_half : temp - p_half;
269  temp = temp / p;
270  return temp * p;
271 }
272 
274 double Round2_to_digit(const double x, const int32_t y) {
275  double exp = pow(10, y);
276  return round(x * exp) / exp;
277 }
278 
280 double round_to_digit(const double x, const int32_t y) {
281  double exp = pow(10, y);
282  return round(x * exp) / exp;
283 }
284 
286 double Sin(const double x) {
287  return sin(x);
288 }
289 
291 double Tan(const double x) {
292  return tan(x);
293 }
294 
296 double Tan__(const float x) {
297  return tanf(x);
298 }
299 
301 double Truncate(const double x, const int32_t y) {
302  double p = pow((double)10L, y);
303  int64_t temp = x * p;
304  return temp / p;
305 }
306 
308 float Truncate__(const float x, const int32_t y) {
309  float p = powf((float)10L, y);
310  int64_t temp = x * p;
311  return temp / p;
312 }
313 
315 int16_t Truncate__1(const int16_t x, const int32_t y) {
316  if (y >= 0) {
317  return x;
318  }
319  int32_t p = pow((float)10L, std::abs(y));
320  int64_t temp = x / p;
321 #if defined(__powerpc__) && !defined(__CUDACC__)
322  int16_t xx = x;
323  xx += 1;
324  temp = xx;
325  temp -= 1;
326  temp /= p;
327 #endif
328  return temp * p;
329 }
330 
332 int32_t Truncate__2(const int32_t x, const int32_t y) {
333  if (y >= 0) {
334  return x;
335  }
336  int32_t p = pow((float)10L, std::abs(y));
337  int64_t temp = x / p;
338  return temp * p;
339 }
340 
342 int64_t Truncate__3(const int64_t x, const int32_t y) {
343  if (y >= 0) {
344  return x;
345  }
346  int64_t p = pow((double)10L, std::abs(y));
347  int64_t temp = x / p;
348  return temp * p;
349 }
350 
352 bool isNan(const double x) {
353  return std::isnan(x);
354 }
355 
357 bool isNan__(const float x) {
358  return std::isnan(x);
359 }
360 
362 double conv_4326_900913_x(const double x) {
363  return x * 111319.490778;
364 }
365 
367 double conv_4326_900913_y(const double y) {
368  return 6378136.99911 * log(tan(.00872664626 * y + .785398163397));
369 }
370 
393 double distance_in_meters(const double fromlon,
394  const double fromlat,
395  const double tolon,
396  const double tolat) {
397  double latitudeArc = (fromlat - tolat) * 0.017453292519943295769236907684886;
398  double longitudeArc = (fromlon - tolon) * 0.017453292519943295769236907684886;
399  double latitudeH = sin(latitudeArc * 0.5);
400  latitudeH *= latitudeH;
401  double lontitudeH = sin(longitudeArc * 0.5);
402  lontitudeH *= lontitudeH;
403  double tmp = cos(fromlat * 0.017453292519943295769236907684886) *
404  cos(tolat * 0.017453292519943295769236907684886);
405  return 6372797.560856 * (2.0 * asin(sqrt(latitudeH + tmp * lontitudeH)));
406 }
407 
409 double distance_in_meters__(const float fromlon,
410  const float fromlat,
411  const float tolon,
412  const float tolat) {
413  float latitudeArc = (fromlat - tolat) * 0.017453292519943295769236907684886;
414  float longitudeArc = (fromlon - tolon) * 0.017453292519943295769236907684886;
415  float latitudeH = sinf(latitudeArc * 0.5);
416  latitudeH *= latitudeH;
417  float lontitudeH = sinf(longitudeArc * 0.5);
418  lontitudeH *= lontitudeH;
419  float tmp = cosf(fromlat * 0.017453292519943295769236907684886) *
420  cosf(tolat * 0.017453292519943295769236907684886);
421  return 6372797.560856 * (2.0 * asinf(sqrtf(latitudeH + tmp * lontitudeH)));
422 }
423 
425 double approx_distance_in_meters(const float fromlon,
426  const float fromlat,
427  const float tolon,
428  const float tolat) {
429 #ifdef __CUDACC__
430  float latitudeArc = (fromlat - tolat) * 0.017453292519943295769236907684886;
431  float longitudeArc = (fromlon - tolon) * 0.017453292519943295769236907684886;
432  float latitudeH = __sinf(latitudeArc * 0.5);
433  latitudeH *= latitudeH;
434  float lontitudeH = __sinf(longitudeArc * 0.5);
435  lontitudeH *= lontitudeH;
436  float tmp = __cosf(fromlat * 0.017453292519943295769236907684886) *
437  __cosf(tolat * 0.017453292519943295769236907684886);
438  return 6372797.560856 * (2.0 * asinf(__fsqrt_rd(latitudeH + tmp * lontitudeH)));
439 #else
440  return distance_in_meters__(fromlon, fromlat, tolon, tolat);
441 #endif
442 }
443 
445 float rect_pixel_bin(const double val,
446  const double min,
447  const double max,
448  const int32_t numbins,
449  const int32_t dimensionsize) {
451  float numbinsf = float(numbins);
452  return float(int32_t(float((val - min) / (max - min)) * numbinsf)) *
453  float(dimensionsize) / numbinsf;
454 }
455 
457 float rect_pixel_bin_x(const double valx,
458  const double minx,
459  const double maxx,
460  const double rectwidth,
461  const double offsetx,
462  const int32_t imgwidth) {
463  const float imgwidthf = float(imgwidth);
464  const float rectwidthf = float(rectwidth);
465  double min = minx;
466  float offset = offsetx;
467  if (offset != 0) {
468  offset = fmodf(offset, rectwidthf);
469  if (offset > 0) {
470  offset -= rectwidthf;
471  }
472  min += offset * (maxx - minx) / imgwidthf;
473  }
474  return float(int32_t(float((valx - min) / (maxx - min)) * (imgwidthf - offset) /
475  rectwidthf)) *
476  rectwidthf +
477  offset + rectwidthf / 2.0f;
478 }
479 
481 float rect_pixel_bin_y(const double valy,
482  const double miny,
483  const double maxy,
484  const double rectheight,
485  const double offsety,
486  const int32_t imgheight) {
487  const float imgheightf = float(imgheight);
488  const float rectheightf = rectheight;
489  double min = miny;
490  float offset = offsety;
491  if (offset != 0) {
492  offset = fmodf(offset, rectheightf);
493  if (offset > 0) {
494  offset -= rectheightf;
495  }
496  min += offset * (maxy - miny) / imgheightf;
497  }
498  return float(int32_t(float((valy - min) / (maxy - min)) * (imgheightf - offset) /
499  rectheightf)) *
500  rectheightf +
501  offset + rectheightf / 2.0f;
502 }
503 
505 float reg_hex_horiz_pixel_bin_x(const double valx,
506  const double minx,
507  const double maxx,
508  const double valy,
509  const double miny,
510  const double maxy,
511  const double hexwidth,
512  const double hexheight,
513  const double offsetx,
514  const double offsety,
515  const int32_t imgwidth,
516  const int32_t imgheight) {
517  const float sqrt3 = 1.7320508075688772;
518  const float imgwidthf = float(imgwidth);
519  const float imgheightf = float(imgheight);
520  const float hexwidthf = float(hexwidth);
521  const float hexheightf = float(hexheight);
522 
523  // expand the bounds of the data according
524  // to the input offsets. This is done because
525  // we also expand the image size according to the
526  // offsets because this algorithm layers the hexagon
527  // bins starting at the bottom left corner
528  double xmin = minx;
529  float xoffset = offsetx;
530  if (xoffset != 0) {
531  xoffset = fmodf(xoffset, hexwidthf);
532  if (xoffset > 0) {
533  xoffset -= hexwidthf;
534  }
535  xmin += xoffset * (maxx - xmin) / imgwidthf;
536  }
537 
538  double ymin = miny;
539  float yoffset = offsety;
540  if (yoffset != 0) {
541  yoffset = fmodf(yoffset, 1.5f * hexheightf);
542  if (yoffset > 0) {
543  yoffset -= 1.5f * hexheightf;
544  }
545  ymin += yoffset * (maxy - ymin) / imgheightf;
546  }
547 
548  // get the pixel position of the point
549  // assumes a linear scale here
550  // Rounds to the nearest pixel.
551  const float pix_x =
552  roundf((imgwidthf - xoffset) * float((valx - xmin) / (maxx - xmin)));
553  const float pix_y =
554  roundf((imgheightf - yoffset) * float((valy - ymin) / (maxy - ymin)));
555 
556  // Now convert the pixel position into a
557  // cube-coordinate system representation
558  const float hexsize = hexheightf / 2.0f;
559  const float cube_x = ((pix_x / sqrt3) - (pix_y / 3.0)) / hexsize;
560  const float cube_z = (pix_y * 2.0f / 3.0f) / hexsize;
561  const float cube_y = -cube_x - cube_z;
562 
563  // need to round the cube coordinates above
564  float rx = round(cube_x);
565  float ry = round(cube_y);
566  float rz = round(cube_z);
567  const float x_diff = fabs(rx - cube_x);
568  const float y_diff = fabs(ry - cube_y);
569  const float z_diff = fabs(rz - cube_z);
570  if (x_diff > y_diff && x_diff > z_diff) {
571  rx = -ry - rz;
572  } else if (y_diff > z_diff) {
573  ry = -rx - rz;
574  } else {
575  rz = -rx - ry;
576  }
577 
578  // now convert the cube/hex coord to a pixel location
579  return hexsize * sqrt3 * (rx + rz / 2.0f) + xoffset;
580 }
581 
583 float reg_hex_horiz_pixel_bin_y(const double valx,
584  const double minx,
585  const double maxx,
586  const double valy,
587  const double miny,
588  const double maxy,
589  const double hexwidth,
590  const double hexheight,
591  const double offsetx,
592  const double offsety,
593  const int32_t imgwidth,
594  const int32_t imgheight) {
595  const float sqrt3 = 1.7320508075688772;
596  const float imgwidthf = float(imgwidth);
597  const float imgheightf = float(imgheight);
598  const float hexwidthf = float(hexwidth);
599  const float hexheightf = float(hexheight);
600 
601  // expand the bounds of the data according
602  // to the input offsets. This is done because
603  // we also expand the image size according to the
604  // offsets because this algorithm layers the hexagon
605  // bins starting at the bottom left corner
606  double xmin = minx;
607  float xoffset = offsetx;
608  if (xoffset != 0) {
609  xoffset = fmodf(xoffset, hexwidthf);
610  if (xoffset > 0) {
611  xoffset -= hexwidthf;
612  }
613  xmin += xoffset * (maxx - xmin) / imgwidthf;
614  }
615 
616  double ymin = miny;
617  float yoffset = offsety;
618  if (yoffset != 0) {
619  yoffset = fmodf(yoffset, 1.5f * hexheightf);
620  if (yoffset > 0) {
621  yoffset -= 1.5f * hexheightf;
622  }
623  ymin += yoffset * (maxy - ymin) / imgheightf;
624  }
625 
626  // get the pixel position of the point
627  // assumes a linear scale here
628  // Rounds to the nearest pixel.
629  const float pix_x =
630  roundf((imgwidthf - xoffset) * float((valx - xmin) / (maxx - xmin)));
631  const float pix_y =
632  roundf((imgheightf - yoffset) * float((valy - ymin) / (maxy - ymin)));
633 
634  // Now convert the pixel position into a
635  // cube-coordinate system representation
636  const float hexsize = hexheightf / 2.0f;
637  const float cube_x = ((pix_x / sqrt3) - (pix_y / 3.0f)) / hexsize;
638  const float cube_z = (pix_y * 2.0f / 3.0f) / hexsize;
639  const float cube_y = -cube_x - cube_z;
640 
641  // need to round the cube coordinates above
642  float rx = round(cube_x);
643  float ry = round(cube_y);
644  float rz = round(cube_z);
645  const float x_diff = fabs(rx - cube_x);
646  const float y_diff = fabs(ry - cube_y);
647  const float z_diff = fabs(rz - cube_z);
648  if ((x_diff <= y_diff || x_diff <= z_diff) && y_diff <= z_diff) {
649  rz = -rx - ry;
650  }
651 
652  // now convert the cube/hex coord to a pixel location
653  return hexsize * 3.0f / 2.0f * rz + yoffset;
654 }
655 
657 float reg_hex_vert_pixel_bin_x(const double valx,
658  const double minx,
659  const double maxx,
660  const double valy,
661  const double miny,
662  const double maxy,
663  const double hexwidth,
664  const double hexheight,
665  const double offsetx,
666  const double offsety,
667  const int32_t imgwidth,
668  const int32_t imgheight) {
669  const float sqrt3 = 1.7320508075688772;
670  const float imgwidthf = float(imgwidth);
671  const float imgheightf = float(imgheight);
672  const float hexwidthf = float(hexwidth);
673  const float hexheightf = float(hexheight);
674 
675  // expand the bounds of the data according
676  // to the input offsets. This is done because
677  // we also expand the image size according to the
678  // offsets because this algorithm layers the hexagon
679  // bins starting at the bottom left corner
680  double xmin = minx;
681  float xoffset = offsetx;
682  if (xoffset != 0) {
683  xoffset = fmodf(xoffset, 1.5f * hexwidthf);
684  if (xoffset > 0) {
685  xoffset -= 1.5f * hexwidthf;
686  }
687  xmin += xoffset * (maxx - xmin) / imgwidthf;
688  }
689 
690  double ymin = miny;
691  float yoffset = offsety;
692  if (yoffset != 0) {
693  yoffset = fmodf(yoffset, hexheightf);
694  if (yoffset > 0) {
695  yoffset -= hexheightf;
696  }
697  ymin += yoffset * (maxy - ymin) / imgheightf;
698  }
699 
700  // get the pixel position of the point
701  // assumes a linear scale here
702  // Rounds to the nearest pixel.
703  const float pix_x =
704  roundf((imgwidthf - xoffset) * float((valx - xmin) / (maxx - xmin)));
705  const float pix_y =
706  roundf((imgheightf - yoffset) * float((valy - ymin) / (maxy - ymin)));
707 
708  // Now convert the pixel position into a
709  // cube-coordinate system representation
710  const float hexsize = hexwidthf / 2.0f;
711  const float cube_x = (pix_x * 2.0f / 3.0f) / hexsize;
712  const float cube_z = ((pix_y / sqrt3) - (pix_x / 3.0f)) / hexsize;
713  const float cube_y = -cube_x - cube_z;
714 
715  // need to round the cube coordinates above
716  float rx = round(cube_x);
717  float ry = round(cube_y);
718  float rz = round(cube_z);
719  const float x_diff = fabs(rx - cube_x);
720  const float y_diff = fabs(ry - cube_y);
721  const float z_diff = fabs(rz - cube_z);
722  if (x_diff > y_diff && x_diff > z_diff) {
723  rx = -ry - rz;
724  }
725 
726  // now convert the cube/hex coord to a pixel location
727  return hexsize * 3.0f / 2.0f * rx + xoffset;
728 }
729 
731 float reg_hex_vert_pixel_bin_y(const double valx,
732  const double minx,
733  const double maxx,
734  const double valy,
735  const double miny,
736  const double maxy,
737  const double hexwidth,
738  const double hexheight,
739  const double offsetx,
740  const double offsety,
741  const int32_t imgwidth,
742  const int32_t imgheight) {
743  const float sqrt3 = 1.7320508075688772;
744  const float imgwidthf = float(imgwidth);
745  const float imgheightf = float(imgheight);
746  const float hexwidthf = float(hexwidth);
747  const float hexheightf = float(hexheight);
748 
749  // expand the bounds of the data according
750  // to the input offsets. This is done because
751  // we also expand the image size according to the
752  // offsets because this algorithm layers the hexagon
753  // bins starting at the bottom left corner
754  float xmin = minx;
755  float xoffset = offsetx;
756  if (xoffset != 0) {
757  xoffset = fmodf(xoffset, 1.5f * hexwidthf);
758  if (xoffset > 0) {
759  xoffset -= 1.5f * hexwidthf;
760  }
761  xmin += xoffset * (maxx - xmin) / imgwidthf;
762  }
763 
764  float ymin = miny;
765  float yoffset = offsety;
766  if (yoffset != 0) {
767  yoffset = fmodf(yoffset, hexheightf);
768  if (yoffset > 0) {
769  yoffset -= hexheightf;
770  }
771  ymin += yoffset * (maxy - ymin) / imgheightf;
772  }
773 
774  // get the pixel position of the point
775  // assumes a linear scale here
776  // Rounds to the nearest pixel.
777  const float pix_x = roundf((imgwidthf - xoffset) * (valx - xmin) / (maxx - xmin));
778  const float pix_y = roundf((imgheightf - yoffset) * (valy - ymin) / (maxy - ymin));
779 
780  // Now convert the pixel position into a
781  // cube-coordinate system representation
782  const float hexsize = hexwidthf / 2.0f;
783  const float cube_x = (pix_x * 2.0f / 3.0f) / hexsize;
784  const float cube_z = ((pix_y / sqrt3) - (pix_x / 3.0f)) / hexsize;
785  const float cube_y = -cube_x - cube_z;
786 
787  // need to round the cube coordinates above
788  float rx = round(cube_x);
789  float ry = round(cube_y);
790  float rz = round(cube_z);
791  const float x_diff = fabs(rx - cube_x);
792  const float y_diff = fabs(ry - cube_y);
793  const float z_diff = fabs(rz - cube_z);
794  if (x_diff > y_diff && x_diff > z_diff) {
795  rx = -ry - rz;
796  } else if (y_diff > z_diff) {
797  ry = -rx - rz;
798  } else {
799  rz = -rx - ry;
800  }
801 
802  // now convert the cube/hex coord to a pixel location
803  return hexsize * sqrt3 * (rz + rx / 2.0f) + yoffset;
804 }
805 
807 double convert_meters_to_merc_pixel_width(const double meters,
808  const double lon,
809  const double lat,
810  const double min_lon,
811  const double max_lon,
812  const int32_t img_width,
813  const double min_width) {
814  const double const1 = 0.017453292519943295769236907684886;
815  const double const2 = 6372797.560856;
816  double t1 = sinf(meters / (2.0 * const2));
817  double t2 = cosf(const1 * lat);
818  const double newlon = lon - (2.0 * asinf(t1 / t2)) / const1;
819  t1 = conv_4326_900913_x(lon);
820  t2 = conv_4326_900913_x(newlon);
821  const double min_merc_x = conv_4326_900913_x(min_lon);
822  const double max_merc_x = conv_4326_900913_x(max_lon);
823  const double merc_diff = max_merc_x - min_merc_x;
824  t1 = ((t1 - min_merc_x) / merc_diff) * static_cast<double>(img_width);
825  t2 = ((t2 - min_merc_x) / merc_diff) * static_cast<double>(img_width);
826 
827  // TODO(croot): need to account for edge cases, such as getting close to the poles.
828  const double sz = fabs(t1 - t2);
829  return (sz < min_width ? min_width : sz);
830 }
831 
833 double convert_meters_to_merc_pixel_height(const double meters,
834  const double lon,
835  const double lat,
836  const double min_lat,
837  const double max_lat,
838  const int32_t img_height,
839  const double min_height) {
840  const double const1 = 0.017453292519943295769236907684886;
841  const double const2 = 6372797.560856;
842  const double latdiff = meters / (const1 * const2);
843  const double newlat =
844  (lat < 0) ? lat + latdiff : lat - latdiff; // assumes a lat range of [-90, 90]
845  double t1 = conv_4326_900913_y(lat);
846  double t2 = conv_4326_900913_y(newlat);
847  const double min_merc_y = conv_4326_900913_y(min_lat);
848  const double max_merc_y = conv_4326_900913_y(max_lat);
849  const double merc_diff = max_merc_y - min_merc_y;
850  t1 = ((t1 - min_merc_y) / merc_diff) * static_cast<double>(img_height);
851  t2 = ((t2 - min_merc_y) / merc_diff) * static_cast<double>(img_height);
852 
853  // TODO(croot): need to account for edge cases, such as getting close to the poles.
854  const double sz = fabs(t1 - t2);
855  return (sz < min_height ? min_height : sz);
856 }
857 
859  const double lat,
860  const double min_lon,
861  const double max_lon,
862  const double min_lat,
863  const double max_lat) {
864  return !(lon < min_lon || lon > max_lon || lat < min_lat || lat > max_lat);
865 }
866 
868  const double lat,
869  const double meters,
870  const double min_lon,
871  const double max_lon,
872  const double min_lat,
873  const double max_lat) {
874  const double const1 = 0.017453292519943295769236907684886;
875  const double const2 = 6372797.560856;
876  const double latdiff = meters / (const1 * const2);
877  const double t1 = sinf(meters / (2.0 * const2));
878  const double t2 = cosf(const1 * lat);
879  const double londiff = (2.0 * asinf(t1 / t2)) / const1;
880  return !(lon + londiff < min_lon || lon - londiff > max_lon ||
881  lat + latdiff < min_lat || lat - latdiff > max_lat);
882 }
883 
884 #include "ExtensionFunctionsGeo.hpp"
EXTENSION_NOINLINE double ln(const double x)
EXTENSION_NOINLINE float reg_hex_horiz_pixel_bin_x(const double valx, const double minx, const double maxx, const double valy, const double miny, const double maxy, const double hexwidth, const double hexheight, const double offsetx, const double offsety, const int32_t imgwidth, const int32_t imgheight)
EXTENSION_NOINLINE double Tan__(const float x)
EXTENSION_INLINE bool is_point_in_merc_view(const double lon, const double lat, const double min_lon, const double max_lon, const double min_lat, const double max_lat)
EXTENSION_NOINLINE double Sin(const double x)
#define EXTENSION_INLINE
EXTENSION_NOINLINE int32_t Truncate__2(const int32_t x, const int32_t y)
EXTENSION_NOINLINE double power(const double x, const double y)
EXTENSION_NOINLINE double Cot(const double x)
EXTENSION_NOINLINE int64_t Round__3(const int64_t x, const int32_t y)
EXTENSION_NOINLINE double distance_in_meters__(const float fromlon, const float fromlat, const float tolon, const float tolat)
EXTENSION_NOINLINE int16_t Round__1(const int16_t x, const int32_t y)
EXTENSION_NOINLINE double Asin(const double x)
EXTENSION_NOINLINE double Log10(const double x)
EXTENSION_NOINLINE double convert_meters_to_merc_pixel_width(const double meters, const double lon, const double lat, const double min_lon, const double max_lon, const int32_t img_width, const double min_width)
EXTENSION_NOINLINE float Floor__(float x)
EXTENSION_NOINLINE int64_t Floor__3(int64_t x)
EXTENSION_NOINLINE double Ceil(double x)
EXTENSION_NOINLINE int16_t Ceil__1(int16_t x)
EXTENSION_NOINLINE int64_t Truncate__3(const int64_t x, const int32_t y)
EXTENSION_NOINLINE int16_t Floor__1(int16_t x)
EXTENSION_NOINLINE double Exp(double x)
EXTENSION_NOINLINE bool isNan__(const float x)
EXTENSION_NOINLINE double radians(const double x)
EXTENSION_NOINLINE double Cos(const double x)
EXTENSION_NOINLINE int64_t Round__4(const int64_t x, const int32_t y0, const int32_t scale)
EXTENSION_NOINLINE int32_t Ceil__2(int32_t x)
EXTENSION_NOINLINE float rect_pixel_bin_y(const double valy, const double miny, const double maxy, const double rectheight, const double offsety, const int32_t imgheight)
EXTENSION_NOINLINE int32_t Round__2(const int32_t x, const int32_t y)
EXTENSION_NOINLINE bool isNan(const double x)
EXTENSION_NOINLINE double Floor(double x)
EXTENSION_NOINLINE double Round2_to_digit(const double x, const int32_t y)
EXTENSION_NOINLINE double degrees(double x)
EXTENSION_NOINLINE double round_to_digit(const double x, const int32_t y)
EXTENSION_NOINLINE double convert_meters_to_merc_pixel_height(const double meters, const double lon, const double lat, const double min_lat, const double max_lat, const int32_t img_height, const double min_height)
EXTENSION_NOINLINE double Log__(const float x)
EXTENSION_NOINLINE float reg_hex_horiz_pixel_bin_y(const double valx, const double minx, const double maxx, const double valy, const double miny, const double maxy, const double hexwidth, const double hexheight, const double offsetx, const double offsety, const int32_t imgwidth, const int32_t imgheight)
EXTENSION_NOINLINE double ln__(const float x)
EXTENSION_NOINLINE float rect_pixel_bin(const double val, const double min, const double max, const int32_t numbins, const int32_t dimensionsize)
EXTENSION_NOINLINE float reg_hex_vert_pixel_bin_x(const double valx, const double minx, const double maxx, const double valy, const double miny, const double maxy, const double hexwidth, const double hexheight, const double offsetx, const double offsety, const int32_t imgwidth, const int32_t imgheight)
#define EXTENSION_NOINLINE
EXTENSION_NOINLINE int64_t Ceil__3(int64_t x)
EXTENSION_NOINLINE float Round__(const float x, const int32_t y)
EXTENSION_NOINLINE float Truncate__(const float x, const int32_t y)
EXTENSION_NOINLINE double distance_in_meters(const double fromlon, const double fromlat, const double tolon, const double tolat)
Computes the distance, in meters, between two WGS-84 positions.
EXTENSION_NOINLINE double Round(const double x, const int32_t y)
EXTENSION_NOINLINE float Ceil__(float x)
EXTENSION_NOINLINE float reg_hex_vert_pixel_bin_y(const double valx, const double minx, const double maxx, const double valy, const double miny, const double maxy, const double hexwidth, const double hexheight, const double offsetx, const double offsety, const int32_t imgwidth, const int32_t imgheight)
EXTENSION_NOINLINE int32_t Floor__2(int32_t x)
EXTENSION_NOINLINE double Truncate(const double x, const int32_t y)
EXTENSION_NOINLINE double Log(const double x)
EXTENSION_NOINLINE double approx_distance_in_meters(const float fromlon, const float fromlat, const float tolon, const float tolat)
EXTENSION_NOINLINE float rect_pixel_bin_x(const double valx, const double minx, const double maxx, const double rectwidth, const double offsetx, const int32_t imgwidth)
EXTENSION_NOINLINE double pi()
EXTENSION_NOINLINE double Tan(const double x)
EXTENSION_NOINLINE double Atan2(const double y, const double x)
EXTENSION_NOINLINE double conv_4326_900913_y(const double y)
EXTENSION_NOINLINE double Atan(const double x)
EXTENSION_NOINLINE double conv_4326_900913_x(const double x)
EXTENSION_NOINLINE double Log10__(const float x)
EXTENSION_NOINLINE bool is_point_size_in_merc_view(const double lon, const double lat, const double meters, const double min_lon, const double max_lon, const double min_lat, const double max_lat)
EXTENSION_NOINLINE int16_t Truncate__1(const int16_t x, const int32_t y)
EXTENSION_NOINLINE double Acos(const double x)