[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[gnuastro-commits] master 4979e34: Table: angular-distance is a new colu
From: |
Mohammad Akhlaghi |
Subject: |
[gnuastro-commits] master 4979e34: Table: angular-distance is a new column arithmetic operator |
Date: |
Wed, 25 Sep 2019 14:51:36 -0400 (EDT) |
branch: master
commit 4979e34ab427a00b4f3bcd2f2538b12871ddc217
Author: Mohammad Akhlaghi <address@hidden>
Commit: Mohammad Akhlaghi <address@hidden>
Table: angular-distance is a new column arithmetic operator
`angular-distance' is a new operator to easily find the distance between
points in various table columns, or the distances of all the rows with a
fixed point.
---
NEWS | 5 ++++
bin/table/arithmetic.c | 76 +++++++++++++++++++++++++++++++++++++++++++++++---
bin/table/arithmetic.h | 1 +
doc/gnuastro.texi | 29 +++++++++++++++++++
4 files changed, 107 insertions(+), 4 deletions(-)
diff --git a/NEWS b/NEWS
index 14c5897..d777388 100644
--- a/NEWS
+++ b/NEWS
@@ -22,6 +22,11 @@ See the end of the file for license conditions.
that have a value of 2, 4 and 5 in the `ID' column.
--notequal: Output only rows that have a different value compared to the
values given to this option in the given column.
+ - Column Arithmetic operators:
+ - `angular-distance': a new operator to easily find the angular
+ distance (along a great circle) between points in various table
+ columns, or the distances of all the points in the table rows with a
+ fixed point. See the book for examples and better explanation.
** Removed features
diff --git a/bin/table/arithmetic.c b/bin/table/arithmetic.c
index 5eb073f..279a9a5 100644
--- a/bin/table/arithmetic.c
+++ b/bin/table/arithmetic.c
@@ -106,6 +106,7 @@ arithmetic_operator_name(int operator)
{
case ARITHMETIC_TABLE_OP_WCSTOIMG: out="wcstoimg"; break;
case ARITHMETIC_TABLE_OP_IMGTOWCS: out="imgtowcs"; break;
+ case ARITHMETIC_TABLE_OP_ANGULARDISTANCE: out="angular-distance"; break;
default:
error(EXIT_FAILURE, 0, "%s: a bug! Please contact us at %s to fix "
"the problem. %d is not a recognized operator code", __func__,
@@ -153,9 +154,11 @@ arithmetic_set_operator(struct tableparams *p, char
*string,
if( op==GAL_ARITHMETIC_OP_INVALID )
{
if( !strncmp(string, "wcstoimg", 8))
- { op=ARITHMETIC_TABLE_OP_WCSTOIMG; *num_operands=0; }
+ { op=ARITHMETIC_TABLE_OP_WCSTOIMG; *num_operands=0; }
else if (!strncmp(string, "imgtowcs", 8))
- { op=ARITHMETIC_TABLE_OP_IMGTOWCS; *num_operands=0; }
+ { op=ARITHMETIC_TABLE_OP_IMGTOWCS; *num_operands=0; }
+ else if (!strncmp(string, "angular-distance", 8))
+ { op=ARITHMETIC_TABLE_OP_ANGULARDISTANCE; *num_operands=0; }
else
{ op=GAL_ARITHMETIC_OP_INVALID; *num_operands=GAL_BLANK_INT; }
}
@@ -430,6 +433,68 @@ arithmetic_wcs(struct tableparams *p, gal_data_t **stack,
int operator)
+static void
+arithmetic_angular_dist(struct tableparams *p, gal_data_t **stack, int
operator)
+{
+ size_t i, j;
+ double *o, *a1, *a2, *b1, *b2;
+ gal_data_t *a, *b, *tmp, *out;
+
+ /* Pop the columns for point `b'.*/
+ tmp=arithmetic_stack_pop(stack, operator);
+ tmp=gal_data_copy_to_new_type_free(tmp, GAL_TYPE_FLOAT64);
+ b=arithmetic_stack_pop(stack, operator);
+ b=gal_data_copy_to_new_type_free(b, GAL_TYPE_FLOAT64);
+ b->next=tmp;
+
+ /* Pop the columns for point `a'.*/
+ tmp=arithmetic_stack_pop(stack, operator);
+ tmp=gal_data_copy_to_new_type_free(tmp, GAL_TYPE_FLOAT64);
+ a=arithmetic_stack_pop(stack, operator);
+ a=gal_data_copy_to_new_type_free(a, GAL_TYPE_FLOAT64);
+ a->next=tmp;
+
+ /* Make sure the sizes are consistant: note that one point can have a
+ single coordinate, but we don't know which one. */
+ if(a->size!=a->next->size)
+ error(EXIT_FAILURE, 0, "the sizes of the third and fourth operands "
+ "of the `%s' operator (respectively containing %zu and %zu "
+ "numbers) must be equal", arithmetic_operator_name(operator),
+ a->next->size, a->size);
+ if(b->size!=b->next->size)
+ error(EXIT_FAILURE, 0, "the sizes of the third and fourth operands "
+ "of the `%s' operator (respectively containing %zu and %zu "
+ "numbers) must be equal", arithmetic_operator_name(operator),
+ b->next->size, b->size);
+
+ /* Make the output array based on the largest size. */
+ out=gal_data_alloc(NULL, GAL_TYPE_FLOAT64, 1,
+ (a->size>b->size ? &a->size : &b->size), NULL, 0,
+ p->cp.minmapsize, p->cp.quietmmap, "angular-dist",
+ NULL, "Angular distance between input points");
+
+ /* Measure the distances. */
+ o=out->array;
+ a1=a->array; a2=a->next->array;
+ b1=b->array; b2=b->next->array;
+ if(a->size==1 || b->size==1) /* One of them is a single point. */
+ for(i=0;i<a->size;++i)
+ for(j=0;j<b->size;++j)
+ o[a->size>b->size?i:j] = gal_wcs_angular_distance_deg(a1[i], a2[i],
+ b1[j], b2[j]);
+ else /* Both have the same length. */
+ for(i=0;i<a->size;++i) /* (all were originally from the same table) */
+ o[i] = gal_wcs_angular_distance_deg(a1[i], a2[i], b1[i], b2[i]);
+
+ /* Clean up and put the output dataset onto the stack. */
+ gal_list_data_free(a);
+ gal_list_data_free(b);
+ gal_list_data_add(stack, out);
+}
+
+
+
+
@@ -545,6 +610,10 @@ arithmetic_operator_run(struct tableparams *p, gal_data_t
**stack,
arithmetic_wcs(p, stack, operator);
break;
+ case ARITHMETIC_TABLE_OP_ANGULARDISTANCE:
+ arithmetic_angular_dist(p, stack, operator);
+ break;
+
default:
error(EXIT_FAILURE, 0, "%s: a bug! Please contact us at %s to "
"fix the problem. The operator code %d is not recognized",
@@ -605,7 +674,7 @@ arithmetic_reverse_polish(struct tableparams *p, struct
column_pack *outpack)
/* A small sanity check. */
if(single->size==1 && p->table && single->size!=p->table->size)
error(EXIT_FAILURE, 0, "the arithmetic operation resulted in a "
- "a single value, but other columns have also been requested "
+ "single value, but other columns have also been requested "
"which have more elements/rows");
/* Set `single->next' to NULL so it isn't treated as a list and
@@ -643,7 +712,6 @@ arithmetic_operate(struct tableparams *p)
size_t i;
struct column_pack *outpack;
-
/* From now on, we will be looking for columns from the index in
`colarray', so to keep things clean, we'll set all the `next' elements
to NULL. */
diff --git a/bin/table/arithmetic.h b/bin/table/arithmetic.h
index dc75d61..7695484 100644
--- a/bin/table/arithmetic.h
+++ b/bin/table/arithmetic.h
@@ -37,6 +37,7 @@ enum arithmetic_operators
{
ARITHMETIC_TABLE_OP_WCSTOIMG = GAL_ARITHMETIC_OP_LAST_CODE,
ARITHMETIC_TABLE_OP_IMGTOWCS,
+ ARITHMETIC_TABLE_OP_ANGULARDISTANCE,
};
diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index 48e4e7a..af7d08a 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -8801,6 +8801,35 @@ $ asttable table.fits -cID,RA -cDEC \
@item imgtowcs
Similar to @code{wcstoimg}, except that image/dataset coordinates are
converted to WCS coordinates.
+
+@item angular-distance
+Return the spherical angular distance (along a great circle, in degrees)
between the given two points.
+Note that each point needs two coordinates (in degrees), so this operator
needs four operands.
+The first and second popped operands are considered to belong to one point and
the third and fourth popped operands to the second point.
+
+Each of the input points can be a single coordinate or a full table column
(containing many points.
+In other words, the following commands are all valid:
+
+@example
+$ asttable table.fits \
+ -c"arith RA1 DEC1 RA2 DEC2 angular-distance"
+$ asttable table.fits \
+ -c"arith RA DEC 12.345 6.789 angular-distance"
+$ asttable table.fits \
+ -c"arith 12.345 6.789 RA DEC angular-distance"
+@end example
+
+In the first case we are assuming that @file{table.fits} has the following
four columns @code{RA1}, @code{DEC1}, @code{RA2}, @code{DEC2}.
+The returned column by this operator would be the difference between two
points in each row with coordinates like the following (@code{RA1},
@code{DEC1}) and (@code{RA2}, @code{DEC2}).
+In other words, for each row, the angular distance between different points is
calculated.
+
+In the second and third cases (which are identical), it is assumed that
@file{table.fits} has the two columns @code{RA} and @code{DEC}.
+The returned column by this operator will be the difference of each row with
the fixed point of (12.345, 6.789).
+
+The distance (along a great circle) on a sphere between two points is
calculated with the equation below, where @mymath{r_1}, @mymath{r_2},
@mymath{d_1} and @mymath{d_2} are the right ascensions and declinations of
points 1 and 2.
+
+@dispmath {\cos(d)=\sin(d_1)\sin(d_2)+\cos(d_1)\cos(d_2)\cos(r_1-r_2)}
+
@end table
[Prev in Thread] |
Current Thread |
[Next in Thread] |
- [gnuastro-commits] master 4979e34: Table: angular-distance is a new column arithmetic operator,
Mohammad Akhlaghi <=