H3C HEALPix library for PostgreSQL  (version 1.2)
h3c_util.c
Go to the documentation of this file.
1 /*
2  Copyright (C) 2012 Gilles Landais (CDS)
3 
4  Author: Gilles Landais, Strasbourg astronomical Data Center (CDS)
5  Email: gilles.landais@unistra.fr
6 
7  This file is part of H3C.
8 
9  H3C is free software; you can redistribute it and/or modify
10  it under the terms of the GNU General Public License as published by
11  the Free Software Foundation; either version 2 of the License, or
12  (at your option) any later version.
13 
14  H3C is distributed in the hope that it will be useful,
15  but WITHOUT ANY WARRANTY; without even the implied warranty of
16  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17  GNU General Public License for more details.
18 
19  You should have received a copy of the GNU General Public License
20  along with H3C; if not, write to the Free Software
21  Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
22 */
23 
36 #include <math.h>
37 #include <sys/time.h>
38 #include <time.h>
39 
40 #include "common.h"
41 #include "h3c_util.h"
42 #include "h3c_math.h"
43 #include "h3cpp_healpix.h"
44 
45 /* PostgreSQL stuff */
46 #include <postgres.h>
47 
48 
49 
50 /*****************************************************************************/
56 /*****************************************************************************/
57 void h3c_log(int level, char *fmt, ...)
58 {
59 #if _H3C_MUTE
60  if (level == H3C_ERROR) {
61  va_list l;
62  char str[BUFSIZ];
63 
64  va_start(l, fmt);
65  vsprintf(str, fmt, l);
66 #ifdef _H3C_MAIN
67  fprintf(stderr, "(ERROR)H3C: %s\n", str);
68  exit(1);
69 #else
70  elog(ERROR, str);
71 #endif
72  }
73 #else
74  static struct timeval tv;
75  static struct timezone tz;
76  static int t0 = -1;
77  int t;
78  va_list l;
79  char str[BUFSIZ];
80 
81  va_start(l, fmt);
82  vsprintf(str, fmt, l);
83 
84  if (level < 0 || t0 == -1) {
85 #ifdef _H3C_MAIN
86  fprintf(stderr, "(debug)H3C: begin time stat\n");
87 #else
88  elog(DEBUG1, "H3C: begin time stat");
89 #endif
90  gettimeofday(&tv, &tz);
91  t0 = tv.tv_sec*1000000+tv.tv_usec;
92  if (level < 0) return;
93  }
94 
95  gettimeofday(&tv, &tz);
96  t = (tv.tv_sec*1000000+tv.tv_usec)-t0;
97 #ifdef _H3C_MAIN
98  fprintf(stderr, "(%s)H3C: %s [%.3f msec]\n",
99  (level==H3C_DEBUG)?"debug":(level==H3C_INFO)?"info":"error",
100  str, t/1000.);
101  if (level == H3C_ERROR) exit(1);
102 #else
103  elog(level, "H3C: %s [%.3f msec]",
104  str, t/1000.);
105 #endif
106  va_end(l);
107 #endif
108 }
109 
110 /*****************************************************************************/
113 /*****************************************************************************/
115 {
116  h3c_log(-1, "chrono");
117 }
118 
119 /*****************************************************************************/
124 /*****************************************************************************/
125 h3c_coord_t h3c_theta2dec(h3c_coord_t theta)
126 {
127  return H3C_PI_2-theta;
128 }
129 
130 /*****************************************************************************/
135 /*****************************************************************************/
136 h3c_coord_t h3c_dec2theta(h3c_coord_t dec)
137 {
138  if (dec>H3C_PI_2) dec=H3C_PI_2;
139  else if(dec<-1*H3C_PI_2) dec=-1*H3C_PI_2;
140  return H3C_PI_2-dec;
141 }
142 
143 /*****************************************************************************/
148 /*****************************************************************************/
149 h3c_coord_t h3c_phi2ra(h3c_coord_t phi)
150 {
151  if (phi > H3C_PI) return phi-H3C_2PI;
152  return phi;
153 }
154 
155 /*****************************************************************************/
160 /*****************************************************************************/
161 h3c_coord_t h3c_ra2phi(h3c_coord_t ra)
162 {
163  if (ra < 0) return H3C_2PI+ra;
164  return ra;
165 }
166 
167 /*****************************************************************************/
172 /*****************************************************************************/
173 int h3c_order(int nside) {
174  unsigned int res=0;
175  if (nside <= 0 || (nside & (nside-1))) {
176  h3c_log(H3C_ERROR, "unvalid nside");
177  return -1;
178  }
179 
180  while (nside > 0x0000FFFF) { res+=16; nside>>=16; }
181  if (nside > 0x000000FF) { res|=8; nside>>=8; }
182  if (nside > 0x0000000F) { res|=4; nside>>=4; }
183  if (nside > 0x00000003) { res|=2; nside>>=2; }
184  if (nside > 0x00000001) { res|=1; }
185  return res;
186 }
187 
188 /*****************************************************************************/
193 /*****************************************************************************/
194 h3c_ipix_t h3c_npix(int nside) {
195  return 12*nside*nside;
196 }
197 
198 /*****************************************************************************/
203 /*****************************************************************************/
204 h3c_coord_t h3c_radius(int nside) {
205  _h3c_vector a, b, c;
206  double t;
207  h3c_unit_vector(2./3., H3C_PI/(4*nside), a);
208  t = 1.-(1./nside);
209  t *= t;
210  h3c_unit_vector(1-t/3., 0, b);
211  h3c_vect_prod(a, b, c);
212 
213  return atan2 (h3c_vect_length(c), h3c_scalar_product(a, b))*H3C_RADEG;
214 }
215 
216 /*****************************************************************************/
224 /*****************************************************************************/
225 h3c_ipix_t *h3c_disk_ipix(h3c_coord_t ra,
226  h3c_coord_t dec,
227  h3c_coord_t radius,
228  int nside,
229  int *count)
230 {
231  return h3cpp_disk_ipix(ra, dec, radius, nside, count, 1);
232 }
233 
234 #if 0
235 /* not effiscient because get too much ipix */
236 
237 #define H3C_MAX_MEM_IPIX 6
238 
239 /*****************************************************************************/
240 void h3c_sorted_list_insert(h3c_ipix_t *list, int count, h3c_ipix_t ipix) {
241  int i;
242  list[count] = ipix;
243  for (i = count-1; i > -1; i--) {
244  if (list[i] < ipix) return;
245  list[i+1] = list[i];
246  list[i] = ipix;
247  }
248 }
249 
250 /*****************************************************************************/
251 /* idea: get the ipix corresponing to the coordinates and search recursively
252  * all ipix arround (according to the property:
253  * dist(MO) < radius+h3c_radius(nside) (O=(ra,dec), M=centre of the ipix)
254  */
255 h3c_ipix_t *h3c_disk_ipix_1(h3c_coord_t ra,
256  h3c_coord_t dec,
257  h3c_coord_t radius,
258  int nside,
259  int *count)
260 {
261  h3c_ipix_t ipix0;
262  static h3c_ipix_t ipix[H3C_MAX_MEM_IPIX]; /* should never over H3C_MAX_MEM_IPIX */
263  h3c_ipix_t *near; /* neighbors */
264  h3c_coord_t theta,phi, ra0, dec0;
265  int i;
266 
267  /* get the ipix where is the center of the disk */
268  h3c_ang2pix_nest(nside,
269  h3c_dec2theta(dec*H3C_DEGRA),
270  h3c_ra2phi(ra*H3C_DEGRA),
271  (h3c_ipix_t *) &ipix0);
272  ipix[0] = ipix0;
273  *count = 1;
274 
275  near = h3cpp_get_neighbors(ipix0, nside);
276  for (i = 0; i < 8 ; i++) {
277  h3c_pix2ang_nest(nside,
278  near[i],
279  &theta,
280  &phi);
281  dec0 = h3c_theta2dec(theta)*H3C_RADEG;
282  ra0 = h3c_phi2ra(phi)*H3C_RADEG;
283  if (h3c_dist(ra, dec, ra0, dec0) < radius) {
284  h3c_log(H3C_DEBUG,"ADD ipix %lld", near[i]);
285  if ((*count) > H3C_MAX_MEM_IPIX-1) h3c_log(H3C_ERROR, "too much ipix!");
286  h3c_sorted_list_insert(ipix, *count, near[i]);
287  (*count)++;
288  }
289  }
290 
291  *count = (*count > 4) ? 4: *count; /* some ipix are forgotten !!!*/
292  return ipix;
293 }
294 #endif
295 
296 /*****************************************************************************/
297 /*****************************************************************************/
298 static int h3c_ipix_t_cmp(const void *a, const void *b) {
299  h3c_ipix_t t =(*(h3c_ipix_t*)a)-*((h3c_ipix_t*)b);
300  if (t==0) return 0;
301  if (t>0) return 1;
302  return -1;
303 /* return (*(h3c_ipix_t*)a)-*((h3c_ipix_t*)b);*/
304 }
305 
306 /*****************************************************************************/
310 /*****************************************************************************/
311 static int h3c_split_list(h3c_ipix_t *in,
312  int nin,
313  h3c_ipix_t *out,
314  int nout,
315  int step)
316 {
317  int i, j = 0;
318  out[0] = in[0];
319  for (i = 1; i < nin; i++) {
320  if (in[i]-in[i-1] > step) {
321  j += 2;
322  if (j < nout) {
323  out[j-1] = in[i-1];
324  out[j] = in[i];
325  }
326  else return -1;
327  }
328  }
329  if (j+1 < nout) out[++j] = in[i-1];
330  else return -1;
331  return j+1;
332 }
333 
334 /*****************************************************************************/
349 /*****************************************************************************/
350 #ifndef LIMIT_CIRCLE_ORDER_0
351 #define LIMIT_CIRCLE_ORDER_0 2/60.
352 #endif
353 #ifndef LIMIT_CIRCLE_ORDER_1
354 #define LIMIT_CIRCLE_ORDER_1 0.5
355 #endif
356 #ifndef LIMIT_CIRCLE_ORDER_2
357 #define LIMIT_CIRCLE_ORDER_2 1
358 #endif
359 #ifndef LIMIT_CIRCLE_ORDER_3
360 #define LIMIT_CIRCLE_ORDER_3 15
361 #endif
362 #ifndef LIMIT_CIRCLE_ORDER_4
363 #define LIMIT_CIRCLE_ORDER_4 40
364 #endif
365 
366 int h3c_disk_ipix_it(h3c_coord_t ra,
367  h3c_coord_t dec,
368  h3c_coord_t radius,
369  int nside,
370  h3c_ipix_t *list,
371  int listsize)
372 {
373  int j, count, step, gap, nside0 = nside;
374  h3c_ipix_t *all;
375  h3c_log(H3C_DEBUG, "begin getting ipix %d", nside);
376 
377 #if 0
378  all = h3c_query_disc(ra, dec, radius, nside, &count);
379 #else
380  /* redefine nside in function of radius */
381  if (radius >= LIMIT_CIRCLE_ORDER_0) {
382  if (radius >= LIMIT_CIRCLE_ORDER_4) nside0 = nside>>6;
383  else if (radius >= LIMIT_CIRCLE_ORDER_3) nside0 = nside>>5;
384  else if (radius >= LIMIT_CIRCLE_ORDER_2) nside0 = nside>>4;
385  else if (radius >= LIMIT_CIRCLE_ORDER_1) nside0 = nside>>3;
386  else nside0 = nside>>2;
387 
388  if (nside0 < 512) nside0 = 512;
389  h3c_log(H3C_INFO, "decrease nside to %d", nside0);
390  }
391 
392  all = h3cpp_disk_ipix(ra, dec, radius, nside0, &count, 1);
393 #endif
394 
395  h3c_log(H3C_DEBUG, "end getting ipix");
396  if (!all) {
397  if (count > 0 ) h3c_log(H3C_ERROR, "index H3C uses too much ipix !");
398  list[0] = -1;
399  return 0;
400  }
401 
402  qsort(all, count, sizeof(h3c_ipix_t), h3c_ipix_t_cmp);
403  h3c_log(H3C_DEBUG, "end sorting");
404 
405  if (radius >= 1) {
406  if (radius >= 10) step = 256;
407  else if (radius >= 5) step = 64;
408  else step = 16;
409  }
410  else step = 4;
411  gap = step;
412 
413  while ( (j = h3c_split_list(all, count, list, listsize, gap)) < 0 ) {
414  gap += step;
415  /*gap *= 2;*/
416  }
417  h3c_log(H3C_DEBUG, "end split");
418 
419  if (nside0 != nside) {
420  for (j=0; j<listsize; j++) {
421  if (j%2) list[j] = (list[j]+1)<<((h3c_order(nside)-h3c_order(nside0))*2);
422  else list[j] = (list[j])<<((h3c_order(nside)-h3c_order(nside0))*2);
423  }
424  }
425 
426  free(all);
427 
428 #if _H3C_MUTE
429 #else
430 { int i, e = 0;
431  static int max = 0;
432  static long long sum = 0;
433  static int npass;
434 
435  for (i = 0; i < j/2; i++) e += list[2*i+1]-list[2*i]+1;
436  h3c_log(H3C_DEBUG, "STAT (h3c_disk_ipix_it) #ipix=%d =split=> %d (gap=%d)",
437  count, e, gap);
438 
439  npass++;
440  sum += count;
441  max = (count>max) ? count : max;
442  if (npass % 10000 == 0) h3c_log(H3C_INFO, "npass=%d moy(#ipix)=%f max(#ipix)=%d", npass, (1.*sum/npass), max);
443 }
444 #endif
445  return j;
446 }
447 
448 
449 /*****************************************************************************/
459 /*****************************************************************************/
460 h3c_coord_t h3c_dist(h3c_coord_t ra1,
461  h3c_coord_t dec1,
462  h3c_coord_t ra2,
463  h3c_coord_t dec2)
464 {
465  h3c_coord_t ra1Rad = ra1*H3C_DEGRA;
466  h3c_coord_t ra2Rad = ra2*H3C_DEGRA;
467  h3c_coord_t dec1Rad = dec1*H3C_DEGRA;
468  h3c_coord_t dec2Rad = dec2*H3C_DEGRA;
469 
470  h3c_coord_t d = pow(sin((dec1Rad-dec2Rad)/2.),2);
471  d += pow(sin((ra1Rad-ra2Rad)/2.), 2)*cos(dec1Rad)*cos(dec2Rad);
472  return (2.*asin(sqrt(d)))*H3C_RADEG;
473 }
474 
475 /*****************************************************************************/
485 /*****************************************************************************/
486 h3c_coord_t h3c_sindist(h3c_coord_t ra1,
487  h3c_coord_t dec1,
488  h3c_coord_t ra2,
489  h3c_coord_t dec2)
490 {
491  h3c_coord_t ra1Rad = ra1*H3C_DEGRA;
492  h3c_coord_t ra2Rad = ra2*H3C_DEGRA;
493  h3c_coord_t dec1Rad = dec1*H3C_DEGRA;
494  h3c_coord_t dec2Rad = dec2*H3C_DEGRA;
495 
496  h3c_coord_t d = pow(sinl((dec1Rad-dec2Rad)/2.),2);
497  d += pow(sinl((ra1Rad-ra2Rad)/2.), 2)*cosl(dec1Rad)*cos(dec2Rad);
498  return d ;
499 }
500 
501 /*****************************************************************************/
502 /*****************************************************************************/
503 static h3c_coord_t h3c_polygon_radius(h3c_coord_t *in_ra,
504  h3c_coord_t *in_dec,
505  int n)
506 {
507  int i;
508  h3c_coord_t center[2];
509  h3c_coord_t d, dist = 0;
510 
511  h3c_poly_center(n, in_ra, in_dec, center);
512  for (i=0; i<n; i++) {
513  d = h3c_dist(center[0], center[1], in_ra[i], in_dec[i]);
514  if (d > dist) dist = d;
515  }
516  return dist;
517 }
518 
519 /*****************************************************************************/
531 /*****************************************************************************/
532 #ifndef LIMIT_POLYGONE_ORDER_0
533 #define LIMIT_POLYGONE_ORDER_0 1/60.
534 #endif
535 #ifndef LIMIT_POLYGONE_ORDER_1
536 #define LIMIT_POLYGONE_ORDER_1 0.5
537 #endif
538 #ifndef LIMIT_POLYGONE_ORDER_2
539 #define LIMIT_POLYGONE_ORDER_2 1
540 #endif
541 
542 int h3c_polygon_ipix_it(h3c_coord_t *in_ra,
543  h3c_coord_t *in_dec,
544  int n,
545  int nside,
546  h3c_ipix_t *list,
547  int listsize)
548 {
549  int j, count, gap, step;
550  h3c_ipix_t *all;
551  int nside0 = nside;
552 
553  h3c_coord_t sz = h3c_polygon_radius(in_ra, in_dec, n);
554 h3c_log(H3C_INFO, "size=%f",sz);
555  if (sz >= LIMIT_POLYGONE_ORDER_0) {
556  if (sz >= LIMIT_POLYGONE_ORDER_2) nside0 = nside>>4;
557  else if (sz >= LIMIT_POLYGONE_ORDER_1) nside0 = nside>>3;
558  else nside0 = nside>>2;
559 
560  if (nside0 < 64) nside0 = 64;
561  h3c_log(H3C_INFO, "decrease nside to %d (circumscribed circle r=%g)",
562  nside0, sz);
563  }
564 
565  h3c_log(H3C_DEBUG, "begin getting ipix");
566  all = h3c_polygon_ipix(in_ra, in_dec, n, nside0, &count);
567  if (!all) {
568  list[0] = 0;
569  h3c_log(H3C_INFO,"no healpix found ?..");
570  list[1] = h3c_npix(nside);
571  return 2;
572  }
573 
574  qsort(all, count, sizeof(h3c_ipix_t), h3c_ipix_t_cmp);
575  h3c_log(H3C_DEBUG, "end sorting (%d items)", count);
576 
577  if (count > 1000000) step = 4096;
578  else if (count > 100000) step = 1024;
579  else if (count > 10000) step = 32;
580  else step = 4;
581  gap = step;
582 
583  while ( (j = h3c_split_list(all, count, list, listsize, gap))<0 ) {
584  gap += step;
585  }
586  h3c_log(H3C_DEBUG, "end split");
587 
588  free(all);
589 
590  if (nside0 != nside) {
591  int i;
592  for (i=0; i<listsize/2; i++) {
593  list[2*i+1] = (list[2*i+1]+1)<<((h3c_order(nside)-h3c_order(nside0))*2);
594  list[2*i] = (list[2*i])<<((h3c_order(nside)-h3c_order(nside0))*2);
595  }
596  }
597 
598 #if _H3C_MUTE
599 #else
600 { int i, e = 0;
601  for (i = 0; i < j/2; i++) e += (list[2*i+1]-list[2*i])+1;
602  h3c_log(H3C_DEBUG, "STAT (h3c_polygon_ipix_it) #ipix=%d =split=> %d (gap=%d)",
603  count, e, gap);
604 }
605 #endif
606  return j;
607 }
608 
609 /*****************************************************************************/
614 /*****************************************************************************/
615 char h3c_in_ellipse(h3c_coord_t alpha, h3c_coord_t delta0,
616  h3c_coord_t alpha1, h3c_coord_t delta01, h3c_coord_t d0,
617  h3c_coord_t e, h3c_coord_t PA0)
618 {
619  h3c_coord_t d_alpha = (alpha1 - alpha) * H3C_DEGRA;
620  h3c_coord_t delta1 = delta01 * H3C_DEGRA;
621  h3c_coord_t delta = delta0 * H3C_DEGRA;
622  h3c_coord_t PA = PA0 * H3C_DEGRA;
623  h3c_coord_t d = d0 * H3C_DEGRA;
624 
625  h3c_coord_t t1 = cos(d_alpha);
626  h3c_coord_t t22 = sin(d_alpha);
627  h3c_coord_t t3 = cos(delta1);
628  h3c_coord_t t32 = sin(delta1);
629  h3c_coord_t t6 = cos(delta);
630  h3c_coord_t t26 = sin(delta);
631  h3c_coord_t t9 = cos(d);
632  h3c_coord_t t55 = sin(d);
633 
634  h3c_coord_t t2;
635  h3c_coord_t t4;
636  h3c_coord_t t5;
637  h3c_coord_t t7;
638  h3c_coord_t t8;
639  h3c_coord_t t10;
640  h3c_coord_t t11;
641  h3c_coord_t t13;
642  h3c_coord_t t14;
643  h3c_coord_t t15;
644  h3c_coord_t t18;
645  h3c_coord_t t19;
646  h3c_coord_t t24;
647  h3c_coord_t t31;
648  h3c_coord_t t36;
649  h3c_coord_t t37;
650  h3c_coord_t t38;
651  h3c_coord_t t45;
652 
653  h3c_coord_t t56;
654  h3c_coord_t t57;
655  h3c_coord_t t60;
656  h3c_coord_t t61;
657  h3c_coord_t t63;
658 
659  if ((t3 * t6 * t1 + t32 * t26) < 0)
660  {
661  return 0;
662  }
663 
664  t2 = t1*t1;
665 
666  t4 = t3*t3;
667  t5 = t2*t4;
668 
669  t7 = t6*t6;
670  t8 = t5*t7;
671 
672  t10 = t9*t9;
673  t11 = t7*t10;
674  t13 = cos(PA);
675  t14 = t13*t13;
676  t15 = t14*t10;
677  t18 = t7*t14;
678  t19 = t18*t10;
679 
680  t24 = sin(PA);
681 
682  t31 = t1*t3;
683 
684  t36 = 2.0*t31*t32*t26*t6;
685  t37 = t31*t32;
686  t38 = t26*t6;
687  t45 = t4*t10;
688 
689  t56 = t55*t55;
690  t57 = t4*t7;
691  t60 = -t8+t5*t11+2.0*t5*t15-t5*t19-2.0*t1*t4*t22*t10*t24*t13*t26-t36+2.0*t37*t38*t10-2.0*t37*t38*t15-t45*t14-t45*t2+2.0*t22*t3*t32*t6*t24*t10*t13-t56+t7-t11+t4-t57+t57*t10+t19-t18*t45;
692  t61 = e*e;
693  t63 = t60*t61+t8+t57-t4-t7+t56+t36;
694  return t63 > 0;
695 }
696 
697 /*****************************************************************************/
698 void h3c_pix2ang_nest(long nside,
699  h3c_ipix_t ipnest,
700  h3c_coord_t *theta,
701  h3c_coord_t *phi) {
702  h3cpp_pix2ang_nest(nside, ipnest, theta, phi);
703 }
704 
705 /*****************************************************************************/
706 void h3c_ang2pix_nest(long nside,
707  h3c_coord_t theta,
708  h3c_coord_t phi,
709  h3c_ipix_t *ipnest) {
710  h3cpp_ang2pix_nest(nside, theta, phi, ipnest);
711 }
712 
713 /*****************************************************************************/
714 void h3c_ring2nest(long nside,
715  h3c_ipix_t ipring,
716  h3c_ipix_t *ipnest) {
717  h3cpp_ring2nest(nside, ipring, ipnest);
718 }
719 
720 
h3c_radius
h3c_coord_t h3c_radius(int nside)
get the max radius of ipix for a given nside
Definition: h3c_util.c:204
h3c_sindist
h3c_coord_t h3c_sindist(h3c_coord_t ra1, h3c_coord_t dec1, h3c_coord_t ra2, h3c_coord_t dec2)
calculate the sinus distance (to be compatible with Q3C)
Definition: h3c_util.c:486
h3c_log
void h3c_log(int level, char *fmt,...)
print statistics
Definition: h3c_util.c:57
h3c_dist
h3c_coord_t h3c_dist(h3c_coord_t ra1, h3c_coord_t dec1, h3c_coord_t ra2, h3c_coord_t dec2)
calculate the distance between 2 points
Definition: h3c_util.c:460
h3c_ra2phi
h3c_coord_t h3c_ra2phi(h3c_coord_t ra)
get the phi value from the right ascension
Definition: h3c_util.c:161
h3c_dec2theta
h3c_coord_t h3c_dec2theta(h3c_coord_t dec)
get the theta value from the declination
Definition: h3c_util.c:136
h3c_vect_length
double h3c_vect_length(_h3c_vector v)
get the length of a vector
Definition: h3c_math.c:114
h3c_theta2dec
h3c_coord_t h3c_theta2dec(h3c_coord_t theta)
get the declination from the theta value
Definition: h3c_util.c:125
h3c_debug_init_clock
void h3c_debug_init_clock()
initalize the debugger timer
Definition: h3c_util.c:114
h3cpp_disk_ipix
h3c_ipix_t * h3cpp_disk_ipix(h3c_coord_t ra, h3c_coord_t dec, double radius, int nside, int *nresult, int inclusive)
retrieve the ipix list for a cone
Definition: h3cpp_healpix.cc:49
h3c_phi2ra
h3c_coord_t h3c_phi2ra(h3c_coord_t phi)
get the right ascension from the phi value
Definition: h3c_util.c:149
h3cpp_ring2nest
void h3cpp_ring2nest(long nside, h3c_ipix_t ipring, h3c_ipix_t *ipnest)
transform ipix (RING) to ipix (NEST)
Definition: h3cpp_healpix.cc:220
h3c_npix
h3c_ipix_t h3c_npix(int nside)
get the number of ipix in the entire sphere
Definition: h3c_util.c:194
h3c_disk_ipix
h3c_ipix_t * h3c_disk_ipix(h3c_coord_t ra, h3c_coord_t dec, h3c_coord_t radius, int nside, int *count)
(c++ interface) get the ipix list of a cone
Definition: h3c_util.c:225
h3cpp_get_neighbors
h3c_ipix_t * h3cpp_get_neighbors(h3c_ipix_t ipix, int nside)
get ipix neighbors
Definition: h3cpp_healpix.cc:148
h3cpp_healpix.h
Interface the Healpix C++ library.
h3c_order
int h3c_order(int nside)
get the order from the nside number
Definition: h3c_util.c:173
LIMIT_POLYGONE_ORDER_0
#define LIMIT_POLYGONE_ORDER_0
get the ipix list of a polygone and group them
Definition: h3c_util.c:533
LIMIT_CIRCLE_ORDER_0
#define LIMIT_CIRCLE_ORDER_0
get the ipix list of a cone and group them
Definition: h3c_util.c:351
common.h
h3c_util.h
usefull functions in the library heaplix for PostgreSQL
h3c_unit_vector
void h3c_unit_vector(double z, double phi, _h3c_vector res)
create a unit vector from a z coordinate and an azimuthal angle.
Definition: h3c_math.c:125
h3c_scalar_product
double h3c_scalar_product(_h3c_vector v1, _h3c_vector v2)
scalar product
Definition: h3c_math.c:92
h3c_in_ellipse
char h3c_in_ellipse(h3c_coord_t alpha, h3c_coord_t delta0, h3c_coord_t alpha1, h3c_coord_t delta01, h3c_coord_t d0, h3c_coord_t e, h3c_coord_t PA0)
verify if a point is inside an ellipse
Definition: h3c_util.c:615
h3c_poly_center
void h3c_poly_center(int, h3c_coord_t *, h3c_coord_t *, h3c_coord_t *)
get the polgon center
Definition: h3c_poly_more.c:86
h3c_polygon_ipix
h3c_ipix_t * h3c_polygon_ipix(h3c_coord_t *in_ra, h3c_coord_t *in_dec, int n, int nside, int *count)
retrieve the ipix list for a polygon
Definition: h3c_poly_more.c:426
h3c_vect_prod
void h3c_vect_prod(_h3c_vector v1, _h3c_vector v2, _h3c_vector r)
vectorial product
Definition: h3c_math.c:79
h3c_math.h
mathematique functions