gnuastro-commits
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

[gnuastro-commits] master eef4e80: Arithmetic's binary operators don't n


From: Mohammad Akhlaghi
Subject: [gnuastro-commits] master eef4e80: Arithmetic's binary operators don't need compiled types
Date: Wed, 13 Dec 2017 23:12:02 -0500 (EST)

branch: master
commit eef4e80c315098145a3ef0dda471c6e6daacae4c
Author: Mohammad Akhlaghi <address@hidden>
Commit: Mohammad Akhlaghi <address@hidden>

    Arithmetic's binary operators don't need compiled types
    
    Until now, Arithmetic's binary operators were managed by two source files
    (`lib/arithmetic-binary.c' and `lib/arithmetic-onlyint.c'). Each of these
    managed many operators and when this got added to the number of types, it
    took extremely long to compile Gnuastro.
    
    With this commit, these two files have been removed and for each operator,
    we now have a separate C source file (for example `lib/arithmetic-plus.c'
    or `lib/arithmetic-bitand.c'). This C file only has an almost single line
    function that calls macros defined in
    `lib/gnuastro-internal/arithmetic-binary.h'. There are still many types to
    manage, but since each is complied separately, when compilation is done in
    parallel, this allow us to compile native binary operators for all types in
    a reasonable amount of time.
    
    Therfore, with this commit, there is no more need for the
    `--enable-bin-op-*' configure time options.
---
 NEWS                                               |   7 +
 configure.ac                                       | 188 -----
 doc/gnuastro.texi                                  | 113 +--
 lib/Makefile.am                                    |  45 +-
 .../arithmetic-onlyint.h => arithmetic-and.c}      |  34 +-
 lib/arithmetic-binary.c                            | 491 -----------
 lib/arithmetic-bitand.c                            |  58 ++
 lib/arithmetic-bitlsh.c                            |  58 ++
 lib/arithmetic-bitor.c                             |  58 ++
 lib/arithmetic-bitrsh.c                            |  58 ++
 lib/arithmetic-bitxor.c                            |  58 ++
 lib/arithmetic-divide.c                            |  53 ++
 .../arithmetic-onlyint.h => arithmetic-eq.c}       |  34 +-
 .../arithmetic-onlyint.h => arithmetic-ge.c}       |  34 +-
 .../arithmetic-onlyint.h => arithmetic-gt.c}       |  34 +-
 .../arithmetic-onlyint.h => arithmetic-le.c}       |  34 +-
 .../arithmetic-onlyint.h => arithmetic-lt.c}       |  34 +-
 lib/arithmetic-minus.c                             |  53 ++
 lib/arithmetic-modulo.c                            |  58 ++
 lib/arithmetic-multiply.c                          |  53 ++
 .../arithmetic-onlyint.h => arithmetic-ne.c}       |  34 +-
 lib/arithmetic-onlyint.c                           | 492 -----------
 .../arithmetic-onlyint.h => arithmetic-or.c}       |  34 +-
 lib/arithmetic-plus.c                              |  53 ++
 lib/arithmetic.c                                   | 899 ++++++++++-----------
 .../{arithmetic-onlyint.h => arithmetic-and.h}     |  14 +-
 lib/gnuastro-internal/arithmetic-binary.h          | 250 +++++-
 .../{arithmetic-binary.h => arithmetic-bitand.h}   |  10 +-
 .../{arithmetic-binary.h => arithmetic-bitlsh.h}   |  10 +-
 .../{arithmetic-binary.h => arithmetic-bitor.h}    |  10 +-
 .../{arithmetic-binary.h => arithmetic-bitrsh.h}   |  10 +-
 .../{arithmetic-binary.h => arithmetic-bitxor.h}   |  10 +-
 .../{arithmetic-binary.h => arithmetic-divide.h}   |  10 +-
 .../{arithmetic-binary.h => arithmetic-eq.h}       |  10 +-
 .../{arithmetic-binary.h => arithmetic-ge.h}       |  10 +-
 .../{arithmetic-binary.h => arithmetic-gt.h}       |  10 +-
 lib/gnuastro-internal/arithmetic-internal.h        |   2 +-
 .../{arithmetic-binary.h => arithmetic-le.h}       |  10 +-
 .../{arithmetic-binary.h => arithmetic-lt.h}       |  10 +-
 .../{arithmetic-binary.h => arithmetic-minus.h}    |  10 +-
 .../{arithmetic-binary.h => arithmetic-modulo.h}   |  10 +-
 .../{arithmetic-binary.h => arithmetic-multiply.h} |  10 +-
 .../{arithmetic-binary.h => arithmetic-ne.h}       |  10 +-
 .../{arithmetic-binary.h => arithmetic-or.h}       |  10 +-
 .../{arithmetic-binary.h => arithmetic-plus.h}     |  10 +-
 45 files changed, 1562 insertions(+), 1941 deletions(-)

diff --git a/NEWS b/NEWS
index e674196..9cbbf0a 100644
--- a/NEWS
+++ b/NEWS
@@ -107,6 +107,13 @@ GNU Astronomy Utilities NEWS                          -*- 
outline -*-
 
 ** Removed features
 
+  Installation: The `--enable-bin-op-*' options that were introduced in
+  Gnuastro 0.3 have been removed. By managing the functions in a better
+  manner (a separate source file for each operator), compilation (when done
+  in parallel) takes about the same time as it took with the default (only
+  four) types before. However, Arithmetic's binary operators can now deal
+  with all ten numeric types natively.
+
   MakeCatalog: `--zeropoint' option doesn't have a short option name any
   more. Previously it was `-z' which was confusing because `-x' and `-y'
   were used to refer to image coordinate positions.
diff --git a/configure.ac b/configure.ac
index dc7a7b7..c1e01f0 100644
--- a/configure.ac
+++ b/configure.ac
@@ -371,194 +371,6 @@ AC_DEFINE_UNQUOTED([CONF_SHOWFMT], [" %-20s"],
 
 
 
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-# The native types for binary arithmetic operations, see the manual for a
-# detailed discussion. The initial list of which types to compile can be
-# determined with the `--enable-bin-op-alltypes' option. If they have a
-# value of 1, they will be compiled. It is possible to disable a previously
-# compiled type with with `--disable-bin-op-TYPENAME', or
-# `--enable-bin-op-TYPENAME=no'.
-binop_alltypes=0
-AC_ARG_ENABLE([bin-op-alltypes],
-              [AS_HELP_STRING([--enable-bin-op-alltypes],
-                    [Allow native binary operations for all types.])],
-             [AS_IF([test "x$enable_bin_op_alltypes" != xno],
-                     [binop_alltypes=1], [binop_alltypes=0])], [])
-AS_IF([test "x$binop_alltypes" != x0],
-      [binop_uint8=1
-       binop_int8=1
-       binop_uint16=1
-       binop_int16=1
-       binop_uint32=1
-       binop_int32=1
-       binop_uint64=1
-       binop_int64=1
-       binop_float32=1
-       binop_float64=1],
-      [binop_uint8=1
-       binop_int8=0
-       binop_uint16=0
-       binop_int16=0
-       binop_uint32=0
-       binop_int32=0
-       binop_uint64=1
-       binop_int64=1
-       binop_float32=1
-       binop_float64=1] )
-
-AC_MSG_CHECKING(compilation of 8-bit unsigned int binary operators)
-AC_ARG_ENABLE([bin-op-uint8],
-              [AS_HELP_STRING([--enable-bin-op-uint8],
-                    [Native binary operators on unsigned 8-bit int.])],
-             [AS_IF([test "x$enable_bin_op_uint8" != xno],
-                     [binop_uint8=1], [binop_uint8=0])], [])
-AS_IF([test "x$binop_uint8" != x0], [binoptprint=yes], [binoptprint=no])
-AC_DEFINE_UNQUOTED([GAL_CONFIG_BIN_OP_UINT8], [$binop_uint8],
-                   [Native binary operations on unsigned 8-bit int.])
-AC_SUBST(HAVE_BIN_OP_UINT8, [$binop_uint8])
-AC_MSG_RESULT($binoptprint)
-
-AC_MSG_CHECKING(compilation of 8-bit signed int binary operators)
-AC_ARG_ENABLE([bin-op-int8],
-              [AS_HELP_STRING([--enable-bin-op-int8],
-                    [Native binary operations on int8 data.])],
-             [AS_IF([test "x$enable_bin_op_int8" != xno],
-                     [binop_int8=1], [binop_int8=0])], [])
-AS_IF([test "x$binop_int8" != x0], [binoptprint=yes], [binoptprint=no])
-AC_DEFINE_UNQUOTED([GAL_CONFIG_BIN_OP_INT8], [$binop_int8],
-                   [Native binary operations on int8 data.])
-AC_SUBST(HAVE_BIN_OP_INT8, [$binop_int8])
-AC_MSG_RESULT($binoptprint)
-
-AC_MSG_CHECKING(compilation of 16-bit unsigned int binary operators)
-AC_ARG_ENABLE([bin-op-uint16],
-              [AS_HELP_STRING([--enable-bin-op-uint16],
-                    [Native binary operators on unsigned short data.])],
-             [AS_IF([test "x$enable_bin_op_uint16" != xno],
-                     [binop_uint16=1], [binop_uint16=0])], [])
-AS_IF([test "x$binop_uint16" != x0], [binoptprint=yes], [binoptprint=no])
-AC_DEFINE_UNQUOTED([GAL_CONFIG_BIN_OP_UINT16], [$binop_uint16],
-                   [Native binary operations on unsigned short data.])
-AC_SUBST(HAVE_BIN_OP_UINT16, [$binop_uint16])
-AC_MSG_RESULT($binoptprint)
-
-AC_MSG_CHECKING(compilation of 16-bit signed int binary operators)
-AC_ARG_ENABLE([bin-op-int16],
-              [AS_HELP_STRING([--enable-bin-op-int16],
-                    [Native binary operations on int16 data.])],
-             [AS_IF([test "x$enable_bin_op_int16" != xno],
-                     [binop_int16=1], [binop_int16=0])], [])
-AS_IF([test "x$binop_int16" != x0], [binoptprint=yes], [binoptprint=no])
-AC_DEFINE_UNQUOTED([GAL_CONFIG_BIN_OP_INT16], [$binop_int16],
-                   [Native binary operations on int16 data.])
-AC_SUBST(HAVE_BIN_OP_INT16, [$binop_int16])
-AC_MSG_RESULT($binoptprint)
-
-AC_MSG_CHECKING(compilation of 32-bit unsigned int binary operators)
-AC_ARG_ENABLE([bin-op-uint32],
-              [AS_HELP_STRING([--enable-bin-op-uint32],
-                    [Native binary operators on unsigned int32 data.])],
-             [AS_IF([test "x$enable_bin_op_uint32" != xno],
-                     [binop_uint32=1], [binop_uint32=0])], [])
-AS_IF([test "x$binop_uint32" != x0], [binoptprint=yes], [binoptprint=no])
-AC_DEFINE_UNQUOTED([GAL_CONFIG_BIN_OP_UINT32], [$binop_uint32],
-                   [Native binary operations on unsigned int32 data.])
-AC_SUBST(HAVE_BIN_OP_UINT32, [$binop_uint32])
-AC_MSG_RESULT($binoptprint)
-
-AC_MSG_CHECKING(compilation of 32-bit signed int binary operators)
-AC_ARG_ENABLE([bin-op-int32],
-              [AS_HELP_STRING([--enable-bin-op-int32],
-                    [Native binary operations on int32 data.])],
-             [AS_IF([test "x$enable_bin_op_int32" != xno],
-                     [binop_int32=1], [binop_int32=0])], [])
-AS_IF([test "x$binop_int32" != x0], [binoptprint=yes], [binoptprint=no])
-AC_DEFINE_UNQUOTED([GAL_CONFIG_BIN_OP_INT32], [$binop_int32],
-                   [Native binary operations on int32 data.])
-AC_SUBST(HAVE_BIN_OP_INT32, [$binop_int32])
-AC_MSG_RESULT($binoptprint)
-
-AC_MSG_CHECKING(compilation of 64-bit unsigned int binary operators)
-AC_ARG_ENABLE([bin-op-uint64],
-              [AS_HELP_STRING([--enable-bin-op-uint64],
-                    [Native binary operators on unsigned long data.])],
-             [AS_IF([test "x$enable_bin_op_uint64" != xno],
-                     [binop_uint64=1], [binop_uint64=0])], [])
-AS_IF([test "x$binop_uint64" != x0], [binoptprint=yes], [binoptprint=no])
-AC_DEFINE_UNQUOTED([GAL_CONFIG_BIN_OP_UINT64], [$binop_uint64],
-                   [Native binary operations on unsigned long data.])
-AC_SUBST(HAVE_BIN_OP_UINT64, [$binop_uint64])
-AC_MSG_RESULT($binoptprint)
-
-AC_MSG_CHECKING(compilation of 64-bit signed int binary operators)
-AC_ARG_ENABLE([bin-op-int64],
-              [AS_HELP_STRING([--enable-bin-op-int64],
-                    [Native binary operations on int64 data.])],
-             [AS_IF([test "x$enable_bin_op_int64" != xno],
-                     [binop_int64=1], [binop_int64=0])], [])
-AS_IF([test "x$binop_int64" != x0], [binoptprint=yes], [binoptprint=no])
-AC_DEFINE_UNQUOTED([GAL_CONFIG_BIN_OP_INT64], [$binop_int64],
-                   [Native binary operations on long64 data.])
-AC_SUBST(HAVE_BIN_OP_INT64, [$binop_int64])
-AC_MSG_RESULT($binoptprint)
-
-AC_MSG_CHECKING(compilation of 32-bit floating point binary operators)
-AC_ARG_ENABLE([bin-op-float32],
-              [AS_HELP_STRING([--enable-bin-op-float32],
-                    [Native binary operations on float32 data.])],
-             [AS_IF([test "x$enable_bin_op_float32" != xno],
-                     [binop_float32=1], [binop_float32=0])], [])
-AS_IF([test "x$binop_float32" != x0], [binoptprint=yes], [binoptprint=no])
-AC_DEFINE_UNQUOTED([GAL_CONFIG_BIN_OP_FLOAT32], [$binop_float32],
-                   [Native binary operations on float32 data.])
-AC_SUBST(HAVE_BIN_OP_FLOAT32, [$binop_float32])
-AC_MSG_RESULT($binoptprint)
-
-AC_MSG_CHECKING(compilation of 64-bit floating point binary operators)
-AC_ARG_ENABLE([bin-op-float64],
-              [AS_HELP_STRING([--enable-bin-op-float64],
-                    [Native binary operations on float64 data.])],
-             [AS_IF([test "x$enable_bin_op_float64" != xno],
-                     [binop_float64=1], [binop_float64=0])], [])
-AS_IF([test "x$binop_float64" != x0], [binoptprint=yes], [binoptprint=no])
-AC_DEFINE_UNQUOTED([GAL_CONFIG_BIN_OP_FLOAT64], [$binop_float64],
-                   [Native binary operations on float64 data.])
-AC_SUBST(HAVE_BIN_OP_FLOAT64, [$binop_float64])
-AC_MSG_RESULT($binoptprint)
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
 # Read arguments about which programs to install. After checking if
 # the argument was actually called, remove any value the user might
 # have given by setting them to "yes" if they are not "no". These
diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index 41ecfc1..5fd62c8 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -4518,74 +4518,6 @@ Do not build or install the program named 
@file{progname}. This is
 very similar to the @option{--enable-progname}, but will build and
 install all the other programs except this one.
 
address@hidden --enable-bin-op-uint8
address@hidden --enable-bin-op-int8
address@hidden --enable-bin-op-uint16
address@hidden --enable-bin-op-int16
address@hidden --enable-bin-op-uint32
address@hidden --enable-bin-op-int32
address@hidden --enable-bin-op-uint64
address@hidden --enable-bin-op-int64
address@hidden --enable-bin-op-float32
address@hidden --enable-bin-op-float64
-Enable the binary data-structure operators to work natively on the
-respective type of data (@option{u} stands for unsigned types, see
address@hidden data types}). Some are compiled by default, to disable them
-(or disable any other type), either run @option{enable-bin-op-TYPE=no}, or
-run @option{--disable-bin-op-TYPE}. The final list of enabled/disabled
-types can be inspected in the outputs of @command{./configure} (close to
-the end).
-
-Binary operators, for example @code{+} or @code{>} (greater than), are some
-of the most common operators to the @ref{Arithmetic} program or the
address@hidden function in @ref{Gnuastro library}. To operate
-most efficiently (as fast as possible without using extra memory or CPU
-resources), it is best to rely on the native types of the input data. For
-example, if you want to add an integer array with a floating point array,
-using the native types, means relying the system's internal type conversion
-for each array element, see @ref{Invoking astarithmetic}. If we don't use
-the native conversion, then the integer array has to be converted to the
-same type as the floating point array to do the conversion. This will
-consume memory (to copy the integer array into a new float array) and CPU
-(integer types need much less processing) resources and ultimately slow
-down the running.
-
-There are many binary operators and in order to have them operate natively
-on of each of the above types, the compiler has to prepare for all the
-different combinations of these types. This can greatly slow down the
address@hidden can also greatly increase the file size of the
-library, from a few hundred kilobytes to a few megabytes.} (when you run
address@hidden). For example, with only one type, @command{make} will
-finish in less than a minute, but if you enable all types, it can take
-roughly half an hour. However, the profits of this one-time investment at
-compilation time will be directly felt (more significantly on large
-images/datasets) each time you run Gnuastro programs or libraries, because
-no internal type conversion will be necessary.
-
-If build time is important for you (mainly developers), disabling shared
-libraries and optimizations (as in @ref{Building and debugging}) is the
-first step to take. If you commonly work with very specific data-types, you
-can enable them (and disable the default types that you don't need) with
-these configuration options. Since the outputs of comparison operators are
address@hidden char} (or @code{uint8_t}) type and most astronomical
-datasets are in single precision (32-bit) floating point (@code{float}),
-the recommended minimum enabled types are @code{uint8} and @code{float32}.
-
-GNU/Linux distribution package managers who compile once, for a large
-audience of users who just download the compiled programs and executables,
-are recommended to enable all types to help their users.
-
address@hidden --enable-bin-op-alltypes
-Enable native binary arithmetic operation on all types, see the description
-above for the various types for a full discussion. As discussed there,
-enabling all types can greatly speed up arithmetic operations on any
-arbitrary dataset, but will also slow down the building time of
-Gnuastro. Recall that in practice this only affects the @ref{Arithmetic}
-program and the @code{gal_arithmetic} library binary operators, nothing
-else. This option is strongly recommended when you are building Gnuastro to
-be included in a package manager of a GNU/Linux distribution (or other
-operating systems).
-
 @item --enable-gnulibcheck
 @cindex GNU C library
 @cindex Gnulib: GNU Portability Library
@@ -10337,33 +10269,14 @@ for(i=0;i<100;++i) out[i]=a[i]+b[i];
 
 @noindent
 Relying on the default C type conversion significantly speeds up the
-processing and also requires less RAM (when using very large
-images). However this great advantage comes at the cost of preparing for
-all the combinations of types while building/compiling Gnuastro. With the
-full list of CFITSIO types, compilation can take roughly half an
-hour. However, some types are not too common, therefore Gnuastro comes with
-a set of configure time options letting you enable certain types for native
-compilation. You can see the full list of @option{--enable-bin-op-XXXX}
-options in @ref{Gnuastro configure options}.
-
-When a type isn't enabled for native binary operations, the input data will
-be internally converted to the smallest, larger type that was enabled. This
-can slow down your processing (which is faster for smaller/integer types)
-and consume more RAM (to copy the new type), so if you often deal with data
-of a specific types, it is much better to make the one-time investment at
-compilation time and reap the benefits each time you run
-Gnuastro/Arithmetic. Note all arithmetic operations are done by
address@hidden function in @ref{Gnuastro library}, so the
-choice of native binary operator types will affect any program (within
-Gnuastro or outside of it) that uses this function (including Arithmetic).
+processing and also requires less RAM (when using very large images).
 
 Some operators can only work on integer types (of any length, for example
 bitwise operators) while others only work on floating point types,
 (currently only the @code{pow} operator). In such cases, if the operand
-type(s) are different an error will be printed and internal conversion
-won't occur. Arithmetic also comes with internal type conversion operators
-which you can use to convert the data into the appropriate type, see
address@hidden operators}.
+type(s) are different, an error will be printed. Arithmetic also comes with
+internal type conversion operators which you can use to convert the data
+into the appropriate type, see @ref{Arithmetic operators}.
 
 @cindex Options
 The hyphen (@command{-}) can be used both to specify options (see
@@ -19077,24 +18990,6 @@ of @code{1}, otherwise it will have a value of 
@code{0}. see
 @ref{Implementation of pthread_barrier} for more.
 @end deffn
 
address@hidden Macro GAL_CONFIG_BIN_OP_UINT8
address@hidden Macro GAL_CONFIG_BIN_OP_INT8
address@hidden Macro GAL_CONFIG_BIN_OP_UINT16
address@hidden Macro GAL_CONFIG_BIN_OP_INT16
address@hidden Macro GAL_CONFIG_BIN_OP_UINT32
address@hidden Macro GAL_CONFIG_BIN_OP_INT32
address@hidden Macro GAL_CONFIG_BIN_OP_UINT64
address@hidden Macro GAL_CONFIG_BIN_OP_INT64
address@hidden Macro GAL_CONFIG_BIN_OP_FLOAT32
address@hidden Macro GAL_CONFIG_BIN_OP_FLOAT64
-If binary arithmetic operators were configured for any type, the respective
-macro will have a value of @code{1} (one), otherwise its value will be
address@hidden (zero). Please see the similar configure-time options in
address@hidden configure options} for a thorough explanation. These are only
-relevant for you if you intend to use the binary operators of
address@hidden on datasets}
address@hidden deffn
-
 @cindex 32-bit
 @cindex 64-bit
 @cindex bit-32
diff --git a/lib/Makefile.am b/lib/Makefile.am
index 94319df..0bef380 100644
--- a/lib/Makefile.am
+++ b/lib/Makefile.am
@@ -51,11 +51,15 @@ libgnuastro_la_LIBADD = 
$(top_builddir)/bootstrapped/lib/libgnu.la
 
 
 # Specify the library .c files
-libgnuastro_la_SOURCES = arithmetic.c arithmetic-binary.c                \
-  arithmetic-onlyint.c binary.c blank.c box.c checkset.c convolve.c      \
-  cosmology.c data.c fits.c git.c interpolate.c list.c match.c options.c \
-  permutation.c polygon.c qsort.c dimension.c statistics.c table.c       \
-  tableintern.c threads.c tile.c timing.c txt.c type.c wcs.c
+libgnuastro_la_SOURCES = arithmetic.c arithmetic-and.c arithmetic-bitand.c \
+  arithmetic-bitlsh.c arithmetic-bitor.c arithmetic-bitrsh.c               \
+  arithmetic-bitxor.c arithmetic-divide.c arithmetic-eq.c arithmetic-ge.c  \
+  arithmetic-gt.c arithmetic-le.c arithmetic-lt.c arithmetic-minus.c       \
+  arithmetic-modulo.c arithmetic-multiply.c arithmetic-ne.c                \
+  arithmetic-or.c arithmetic-plus.c binary.c blank.c box.c checkset.c      \
+  convolve.c cosmology.c data.c fits.c git.c interpolate.c list.c match.c  \
+  options.c permutation.c polygon.c qsort.c dimension.c statistics.c       \
+  table.c tableintern.c threads.c tile.c timing.c txt.c type.c wcs.c
 
 
 
@@ -85,12 +89,21 @@ pkginclude_HEADERS = gnuastro/config.h 
$(headersdir)/arithmetic.h         \
 # will not distributed, so we need to explicitly tell Automake to
 # distribute them here.
 internaldir=$(top_srcdir)/lib/gnuastro-internal
-EXTRA_DIST = gnuastro.pc.in $(headersdir)/README $(internaldir)/README    \
-  $(internaldir)/arithmetic-binary.h $(internaldir)/arithmetic-internal.h \
-  $(internaldir)/arithmetic-onlyint.h $(internaldir)/checkset.h           \
-  $(internaldir)/commonopts.h $(internaldir)/config.h.in                  \
-  $(internaldir)/fixedstringmacros.h $(internaldir)/options.h             \
-  $(internaldir)/tableintern.h $(internaldir)/timing.h
+EXTRA_DIST = gnuastro.pc.in $(headersdir)/README $(internaldir)/README  \
+  $(internaldir)/arithmetic-and.h $(internaldir)/arithmetic-binary.h    \
+  $(internaldir)/arithmetic-bitand.h $(internaldir)/arithmetic-bitlsh.h \
+  $(internaldir)/arithmetic-bitor.h $(internaldir)/arithmetic-bitrsh.h  \
+  $(internaldir)/arithmetic-bitxor.h $(internaldir)/arithmetic-divide.h \
+  $(internaldir)/arithmetic-eq.h $(internaldir)/arithmetic-ge.h         \
+  $(internaldir)/arithmetic-gt.h $(internaldir)/arithmetic-internal.h   \
+  $(internaldir)/arithmetic-le.h $(internaldir)/arithmetic-lt.h         \
+  $(internaldir)/arithmetic-minus.h $(internaldir)/arithmetic-modulo.h  \
+  $(internaldir)/arithmetic-multiply.h $(internaldir)/arithmetic-ne.h   \
+  $(internaldir)/arithmetic-or.h $(internaldir)/arithmetic-plus.h       \
+  $(internaldir)/checkset.h $(internaldir)/commonopts.h                 \
+  $(internaldir)/config.h.in $(internaldir)/fixedstringmacros.h         \
+  $(internaldir)/options.h $(internaldir)/tableintern.h                 \
+  $(internaldir)/timing.h
 
 
 
@@ -115,16 +128,6 @@ gnuastro/config.h: Makefile $(internaldir)/config.h.in
               -e 's|@address@hidden|$(HAVE_LIBGIT2)|g'                  \
               -e 's|@address@hidden|$(HAVE_WCSLIB_VERSION)|g'    \
               -e 's|@address@hidden|$(HAVE_PTHREAD_BARRIER)|g'  \
-              -e 's|@address@hidden|$(HAVE_BIN_OP_UINT8)|g'        \
-              -e 's|@address@hidden|$(HAVE_BIN_OP_INT8)|g'          \
-              -e 's|@address@hidden|$(HAVE_BIN_OP_UINT16)|g'      \
-              -e 's|@address@hidden|$(HAVE_BIN_OP_INT16)|g'        \
-              -e 's|@address@hidden|$(HAVE_BIN_OP_UINT32)|g'      \
-              -e 's|@address@hidden|$(HAVE_BIN_OP_INT32)|g'        \
-              -e 's|@address@hidden|$(HAVE_BIN_OP_UINT64)|g'      \
-              -e 's|@address@hidden|$(HAVE_BIN_OP_INT64)|g'        \
-              -e 's|@address@hidden|$(HAVE_BIN_OP_FLOAT32)|g'    \
-              -e 's|@address@hidden|$(HAVE_BIN_OP_FLOAT64)|g'    \
               -e 's|@address@hidden|$(SIZEOF_LONG)|g'                    \
               -e 's|@address@hidden|$(SIZEOF_SIZE_T)|g'                \
               -e 's|@address@hidden|$(RESTRICT_REPLACEMENT)|g'  \
diff --git a/lib/gnuastro-internal/arithmetic-onlyint.h b/lib/arithmetic-and.c
similarity index 50%
copy from lib/gnuastro-internal/arithmetic-onlyint.h
copy to lib/arithmetic-and.c
index 2a8b3bd..aedc84a 100644
--- a/lib/gnuastro-internal/arithmetic-onlyint.h
+++ b/lib/arithmetic-and.c
@@ -5,7 +5,7 @@ This is part of GNU Astronomy Utilities (Gnuastro) package.
 Original author:
      Mohammad Akhlaghi <address@hidden>
 Contributing author(s):
-Copyright (C) 2015, Free Software Foundation, Inc.
+Copyright (C) 2016, Free Software Foundation, Inc.
 
 Gnuastro is free software: you can redistribute it and/or modify it
 under the terms of the GNU General Public License as published by the
@@ -20,16 +20,32 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __GAL_ARITHMETIC_ONLYINT_H__
-#define __GAL_ARITHMETIC_ONLYINT_H__
+#include <config.h>
 
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
 
-gal_data_t *
-arithmetic_onlyint_binary(int operator, int flags, gal_data_t *lo,
-                          gal_data_t *ro);
+#include <gnuastro/blank.h>
 
-gal_data_t *
-arithmetic_onlyint_bitwise_not(int flags, gal_data_t *in);
+#include <gnuastro-internal/arithmetic-binary.h>
+#include <gnuastro-internal/arithmetic-internal.h>
 
 
-#endif
+
+
+
+/* `BINARY_SET_LT' is defined in `arithmetic-binary.h'. As you see there,
+   this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_and(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  int checkblank=gal_arithmetic_binary_checkblank(l, r);
+
+  BINARY_SET_LT( ARITHMETIC_BINARY_OUT_TYPE_UINT8, && );
+}
diff --git a/lib/arithmetic-binary.c b/lib/arithmetic-binary.c
deleted file mode 100644
index a0a69da..0000000
--- a/lib/arithmetic-binary.c
+++ /dev/null
@@ -1,491 +0,0 @@
-/*********************************************************************
-Arithmetic operations on data structures.
-This is part of GNU Astronomy Utilities (Gnuastro) package.
-
-Original author:
-     Mohammad Akhlaghi <address@hidden>
-Contributing author(s):
-Copyright (C) 2016, Free Software Foundation, Inc.
-
-Gnuastro is free software: you can redistribute it and/or modify it
-under the terms of the GNU General Public License as published by the
-Free Software Foundation, either version 3 of the License, or (at your
-option) any later version.
-
-Gnuastro is distributed in the hope that it will be useful, but
-WITHOUT ANY WARRANTY; without even the implied warranty of
-MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
-General Public License for more details.
-
-You should have received a copy of the GNU General Public License
-along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
-**********************************************************************/
-#include <config.h>
-
-#include <stdio.h>
-#include <errno.h>
-#include <error.h>
-#include <stdlib.h>
-
-#include <gnuastro/blank.h>
-#include <gnuastro/arithmetic.h>
-
-#include <gnuastro-internal/arithmetic-internal.h>
-
-
-
-
-
-/************************************************************************/
-/*************            Native type macros            *****************/
-/************************************************************************/
-#if GAL_CONFIG_BIN_OP_UINT8 == 1
-#define BINARY_LT_IS_UINT8                                         \
-  case GAL_TYPE_UINT8: BINARY_LT_SET(uint8_t);          break;
-#define BINARY_LT_SET_RT_IS_UINT8(LT)                              \
-  case GAL_TYPE_UINT8: BINARY_RT_LT_SET(uint8_t, LT);   break;
-#else
-#define BINARY_LT_IS_UINT8
-#define BINARY_LT_SET_RT_IS_UINT8(LT)
-#endif
-
-
-
-
-
-#if GAL_CONFIG_BIN_OP_INT8 == 1
-#define BINARY_LT_IS_INT8                                          \
-  case GAL_TYPE_INT8: BINARY_LT_SET(int8_t);            break;
-#define BINARY_LT_SET_RT_IS_INT8(LT)                               \
-  case GAL_TYPE_INT8: BINARY_RT_LT_SET(int8_t, LT);     break;
-#else
-#define BINARY_LT_IS_INT8
-#define BINARY_LT_SET_RT_IS_INT8(LT)
-#endif
-
-
-
-
-
-#if GAL_CONFIG_BIN_OP_UINT16 == 1
-#define BINARY_LT_IS_UINT16                                        \
-  case GAL_TYPE_UINT16: BINARY_LT_SET(uint16_t);        break;
-#define BINARY_LT_SET_RT_IS_UINT16(LT)                             \
-  case GAL_TYPE_UINT16: BINARY_RT_LT_SET(uint16_t, LT); break;
-#else
-#define BINARY_LT_IS_UINT16
-#define BINARY_LT_SET_RT_IS_UINT16(LT)
-#endif
-
-
-
-
-
-#if GAL_CONFIG_BIN_OP_INT16 == 1
-#define BINARY_LT_IS_INT16                                         \
-  case GAL_TYPE_INT16: BINARY_LT_SET(int16_t);          break;
-#define BINARY_LT_SET_RT_IS_INT16(LT)                              \
-  case GAL_TYPE_INT16: BINARY_RT_LT_SET(int16_t, LT);   break;
-#else
-#define BINARY_LT_IS_INT16
-#define BINARY_LT_SET_RT_IS_INT16(LT)
-#endif
-
-
-
-
-
-#if GAL_CONFIG_BIN_OP_UINT32 == 1
-#define BINARY_LT_IS_UINT32                                        \
-  case GAL_TYPE_UINT32: BINARY_LT_SET(uint32_t);        break;
-#define BINARY_LT_SET_RT_IS_UINT32(LT)                             \
-  case GAL_TYPE_UINT32: BINARY_RT_LT_SET(uint32_t, LT); break;
-#else
-#define BINARY_LT_IS_UINT32
-#define BINARY_LT_SET_RT_IS_UINT32(LT)
-#endif
-
-
-
-
-
-#if GAL_CONFIG_BIN_OP_INT32 == 1
-#define BINARY_LT_IS_INT32                                         \
-  case GAL_TYPE_INT32: BINARY_LT_SET(int32_t);          break;
-#define BINARY_LT_SET_RT_IS_INT32(LT)                              \
-  case GAL_TYPE_INT32: BINARY_RT_LT_SET(int32_t, LT);   break;
-#else
-#define BINARY_LT_IS_INT32
-#define BINARY_LT_SET_RT_IS_INT32(LT)
-#endif
-
-
-
-
-
-#if GAL_CONFIG_BIN_OP_UINT64 == 1
-#define BINARY_LT_IS_UINT64                                        \
-  case GAL_TYPE_UINT64: BINARY_LT_SET(uint64_t);        break;
-#define BINARY_LT_SET_RT_IS_UINT64(LT)                             \
-  case GAL_TYPE_UINT64: BINARY_RT_LT_SET(uint64_t, LT); break;
-#else
-#define BINARY_LT_IS_UINT64
-#define BINARY_LT_SET_RT_IS_UINT64(LT)
-#endif
-
-
-
-
-
-#if GAL_CONFIG_BIN_OP_INT64 == 1
-#define BINARY_LT_IS_INT64                                         \
-  case GAL_TYPE_INT64: BINARY_LT_SET(int64_t);          break;
-#define BINARY_LT_SET_RT_IS_INT64(LT)                              \
-  case GAL_TYPE_INT64: BINARY_RT_LT_SET(int64_t, LT);   break;
-#else
-#define BINARY_LT_IS_INT64
-#define BINARY_LT_SET_RT_IS_INT64(LT)
-#endif
-
-
-
-
-
-#if GAL_CONFIG_BIN_OP_FLOAT32 == 1
-#define BINARY_LT_IS_FLOAT32                                       \
-  case GAL_TYPE_FLOAT32: BINARY_LT_SET(float);          break;
-#define BINARY_LT_SET_RT_IS_FLOAT32(LT)                            \
-  case GAL_TYPE_FLOAT32: BINARY_RT_LT_SET(float, LT);   break;
-#else
-#define BINARY_LT_IS_FLOAT32
-#define BINARY_LT_SET_RT_IS_FLOAT32(LT)
-#endif
-
-
-
-
-
-#if GAL_CONFIG_BIN_OP_FLOAT64 == 1
-#define BINARY_LT_IS_FLOAT64                                       \
-  case GAL_TYPE_FLOAT64: BINARY_LT_SET(double);         break;
-#define BINARY_LT_SET_RT_IS_FLOAT64(LT)                            \
-  case GAL_TYPE_FLOAT64: BINARY_RT_LT_SET(double, LT);  break;
-#else
-#define BINARY_LT_IS_FLOAT64
-#define BINARY_LT_SET_RT_IS_FLOAT64(LT)
-#endif
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-/************************************************************************/
-/*************              High level macros           *****************/
-/************************************************************************/
-/* Final step to be used by all operators and all types. */
-#define BINARY_OP_OT_RT_LT_SET(OP, OT, RT, LT) {                        \
-    LT lb, *la=l->array;                                                \
-    RT rb, *ra=r->array;                                                \
-    OT ob, *oa=o->array, *of=oa + o->size;                              \
-    if(checkblank)                                                      \
-      {                                                                 \
-        gal_blank_write(&lb, l->type);                                  \
-        gal_blank_write(&rb, r->type);                                  \
-        gal_blank_write(&ob, o->type);                                  \
-        do                                                              \
-          {                                                             \
-            if(lb==lb && rb==rb)/* Both are integers.                */ \
-              *oa = (*la!=lb  && *ra!=rb)  ? *la OP *ra : ob ;          \
-            else if(lb==lb)     /* Only left operand is an integer.  */ \
-              *oa = (*la!=lb  && *ra==*ra) ? *la OP *ra : ob;           \
-            else                /* Only right operand is an integer. */ \
-              *oa = (*la==*la && *ra!=rb)  ? *la OP *ra : ob;           \
-            if(l->size>1) ++la;                                         \
-            if(r->size>1) ++ra;                                         \
-          }                                                             \
-        while(++oa<of);                                                 \
-      }                                                                 \
-    else                                                                \
-      {                                                                 \
-        if(l->size==r->size) do *oa = *la++ OP *ra++; while(++oa<of);   \
-        else if(l->size==1)  do *oa = *la   OP *ra++; while(++oa<of);   \
-        else                 do *oa = *la++ OP *ra;   while(++oa<of);   \
-      }                                                                 \
-  }
-
-
-
-
-/* This is for operators like `&&' and `||', where the right operator is
-   not necessarily read (and thus incremented). */
-#define BINARY_OP_INCR_OT_RT_LT_SET(OP, OT, RT, LT) {                   \
-    LT *la=l->array;                                                    \
-    RT *ra=r->array;                                                    \
-    OT *oa=o->array, *of=oa + o->size;                                  \
-    if(l->size==r->size) do {*oa = *la++ OP *ra; ++ra;} while(++oa<of); \
-    else if(l->size==1)  do {*oa = *la   OP *ra; ++ra;} while(++oa<of); \
-    else                 do  *oa = *la++ OP *ra;        while(++oa<of); \
-  }
-
-
-
-
-
-/* For operators whose type may be any of the given inputs. */
-#define BINARY_OP_RT_LT_SET(OP, RT, LT)                            \
-  if(o->type==l->type)                                             \
-    BINARY_OP_OT_RT_LT_SET(OP, LT, RT, LT)                         \
-  else                                                             \
-    BINARY_OP_OT_RT_LT_SET(OP, RT, RT, LT)
-
-
-
-
-
-/* Left and right types set, choose what to do based on operator. */
-#define BINARY_RT_LT_SET(RT, LT)                                    \
-  switch(operator)                                                  \
-    {                                                               \
-    case GAL_ARITHMETIC_OP_PLUS:                                    \
-      BINARY_OP_RT_LT_SET(          +, RT, LT);                     \
-      break;                                                        \
-    case GAL_ARITHMETIC_OP_MINUS:                                   \
-      BINARY_OP_RT_LT_SET(          -, RT, LT);                     \
-      break;                                                        \
-    case GAL_ARITHMETIC_OP_MULTIPLY:                                \
-      BINARY_OP_RT_LT_SET(          *, RT, LT);                     \
-      break;                                                        \
-    case GAL_ARITHMETIC_OP_DIVIDE:                                  \
-      BINARY_OP_RT_LT_SET(          /, RT, LT);                     \
-      break;                                                        \
-    case GAL_ARITHMETIC_OP_LT:                                      \
-      BINARY_OP_OT_RT_LT_SET(       <, uint8_t, RT, LT);            \
-      break;                                                        \
-    case GAL_ARITHMETIC_OP_LE:                                      \
-      BINARY_OP_OT_RT_LT_SET(      <=, uint8_t, RT, LT);            \
-      break;                                                        \
-    case GAL_ARITHMETIC_OP_GT:                                      \
-      BINARY_OP_OT_RT_LT_SET(       >, uint8_t, RT, LT);            \
-      break;                                                        \
-    case GAL_ARITHMETIC_OP_GE:                                      \
-      BINARY_OP_OT_RT_LT_SET(      >=, uint8_t, RT, LT);            \
-      break;                                                        \
-    case GAL_ARITHMETIC_OP_EQ:                                      \
-      BINARY_OP_OT_RT_LT_SET(      ==, uint8_t, RT, LT);            \
-      break;                                                        \
-    case GAL_ARITHMETIC_OP_NE:                                      \
-      BINARY_OP_OT_RT_LT_SET(      !=, uint8_t, RT, LT);            \
-      break;                                                        \
-    case GAL_ARITHMETIC_OP_AND:                                     \
-      BINARY_OP_INCR_OT_RT_LT_SET( &&, uint8_t, RT, LT);            \
-      break;                                                        \
-    case GAL_ARITHMETIC_OP_OR:                                      \
-      BINARY_OP_INCR_OT_RT_LT_SET( ||, uint8_t, RT, LT);            \
-      break;                                                        \
-    default:                                                        \
-      error(EXIT_FAILURE, 0, "%s: operator code %d not recognized", \
-            "BINARY_RT_LT_SET", operator);                          \
-    }
-
-
-
-
-
-
-/* Left operand type set, see what the right operand type is. */
-#define BINARY_LT_SET(LT)                                          \
-  switch(r->type)                                                  \
-    {                                                              \
-      BINARY_LT_SET_RT_IS_UINT8(LT);                               \
-      BINARY_LT_SET_RT_IS_INT8(LT);                                \
-      BINARY_LT_SET_RT_IS_UINT16(LT);                              \
-      BINARY_LT_SET_RT_IS_INT16(LT);                               \
-      BINARY_LT_SET_RT_IS_UINT32(LT);                              \
-      BINARY_LT_SET_RT_IS_INT32(LT);                               \
-      BINARY_LT_SET_RT_IS_UINT64(LT);                              \
-      BINARY_LT_SET_RT_IS_INT64(LT);                               \
-      BINARY_LT_SET_RT_IS_FLOAT32(LT);                             \
-      BINARY_LT_SET_RT_IS_FLOAT64(LT);                             \
-    default:                                                       \
-      error(EXIT_FAILURE, 0, "%s: type code %d not recognized",    \
-            "BINARY_LT_SET", r->type);                             \
-    }
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-/************************************************************************/
-/*************              Top level function          *****************/
-/************************************************************************/
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro)
-{
-  /* Read the variable arguments. `lo' and `ro' keep the original data, in
-     case their type isn't built (based on configure options are configure
-     time). */
-  int checkblank;
-  int32_t otype, final_otype;
-  size_t out_size, minmapsize;
-  gal_data_t *l, *r, *o=NULL, *tmp_o;
-
-
-  /* Simple sanity check on the input sizes */
-  if( !( (flags & GAL_ARITHMETIC_NUMOK) && (lo->size==1 || ro->size==1))
-      && gal_data_dsize_is_different(lo, ro) )
-    error(EXIT_FAILURE, 0, "%s: the non-number inputs to %s don't have the "
-          "same dimension/size", __func__,
-          gal_arithmetic_operator_string(operator));
-
-
-  /* Set the final output type (independent of which types are
-     compiled). These needs to be done before the call to
-     `gal_arithmetic_convert_to_compiled_type', because that function can
-     free the space of the original data structures, thus we will loose the
-     original data structure information. */
-  final_otype=gal_arithmetic_binary_out_type(operator, lo, ro);
-
-
-  /* Make sure the input arrays have one of the compiled types. From this
-     point on, until the cleaning up section of this function, we won't be
-     using the `lo' and `ro' pointers. */
-  l=gal_arithmetic_convert_to_compiled_type(lo, flags);
-  r=gal_arithmetic_convert_to_compiled_type(ro, flags);
-
-
-  /* Set the output type. For the comparison operators, the output type is
-     either 0 or 1, so we will set the output type to `unsigned char' for
-     efficient memory and CPU usage. Since the number of operators without
-     a fixed output type (like the conditionals) is less, by `default' we
-     will set the output type to `unsigned char', and if any of the other
-     operatrs are given, it will be chosen based on the input types.*/
-  otype=gal_arithmetic_binary_out_type(operator, l, r);
-
-
-  /* Set the output sizes. */
-  minmapsize = ( l->minmapsize < r->minmapsize
-                 ? l->minmapsize : r->minmapsize );
-  out_size = l->size > r->size ? l->size : r->size;
-
-
-  /* If we want inplace output, set the output pointer to one input. Note
-     that the output type can be different from both inputs.  */
-  if(flags & GAL_ARITHMETIC_INPLACE)
-    {
-      if     (l->type==otype && out_size==l->size)   o = l;
-      else if(r->type==otype && out_size==r->size)   o = r;
-    }
-
-
-  /* If the output pointer was not set for any reason, allocate it. For
-     `mmapsize', note that since its `size_t', it will always be
-     Positive. The `-1' that is recommended to give when you want the value
-     in RAM is actually the largest possible memory location. So we just
-     have to choose the smaller minmapsize of the two to decide if the
-     output array should be in RAM or not. */
-  if(o==NULL)
-    o = gal_data_alloc(NULL, otype,
-                       l->size>1 ? l->ndim  : r->ndim,
-                       l->size>1 ? l->dsize : r->dsize,
-                       l->size>1 ? l->wcs   : r->wcs,
-                       0, minmapsize, NULL, NULL, NULL );
-
-
-  /* See if we should check for blanks. When both types are floats, blanks
-     don't need to be checked (the floating point standard will do the job
-     for us). It is also not necessary to check blanks in bitwise
-     operators, but bitwise operators have their own macro
-     (`BINARY_OP_INCR_OT_RT_LT_SET') which doesn' use `checkblanks'.*/
-  checkblank = ((((l->type!=GAL_TYPE_FLOAT32    && l->type!=GAL_TYPE_FLOAT64)
-                  || (r->type!=GAL_TYPE_FLOAT32 && r->type!=GAL_TYPE_FLOAT64))
-                 && (gal_blank_present(l, 1) || gal_blank_present(r, 1)))
-                ? 1 : 0 );
-
-
-  /* Start setting the operator and operands. */
-  switch(l->type)
-    {
-      BINARY_LT_IS_UINT8;
-      BINARY_LT_IS_INT8;
-      BINARY_LT_IS_UINT16;
-      BINARY_LT_IS_INT16;
-      BINARY_LT_IS_UINT32;
-      BINARY_LT_IS_INT32;
-      BINARY_LT_IS_UINT64;
-      BINARY_LT_IS_INT64;
-      BINARY_LT_IS_FLOAT32;
-      BINARY_LT_IS_FLOAT64;
-    default:
-      error(EXIT_FAILURE, 0, "%s: type code %d not recognized",
-            __func__, l->type);
-    }
-
-
-  /* Clean up. Note that if the input arrays can be freed, and any of right
-     or left arrays needed conversion, `BINARY_CONVERT_TO_COMPILED_TYPE'
-     has already freed the input arrays, so only `r' and `l' need
-     freeing. Alternatively, when the inputs shouldn't be freed, the only
-     allocated spaces are the `r' and `l' arrays if their types weren't
-     compiled for binary operations, we can tell this from the pointers: if
-     they are different from the original pointers, they were allocated. */
-  if(flags & GAL_ARITHMETIC_FREE)
-    {
-      if     (o==l)       gal_data_free(r);
-      else if(o==r)       gal_data_free(l);
-      else              { gal_data_free(l); gal_data_free(r); }
-    }
-  else
-    {
-      if(l!=lo)           gal_data_free(l);
-      if(r!=ro)           gal_data_free(r);
-    }
-
-
-  /* The type of the output dataset (`o->type') was chosen from `l' and `r'
-     (copies of the orignal operands but in a compiled type, not
-     necessarily the original `lo' and `ro' data structures). So we need to
-     to get the final type based on the original operands and check if the
-     final output needs changing.
-
-     IMPORTANT: This has to be done after (possibly) freeing the left and
-     right operands because this step can change the `o' pointer which
-     they may depend on (when working in-place). */
-  if( o->type != final_otype )
-    {
-      tmp_o=gal_data_copy_to_new_type(o, final_otype);
-      gal_data_free(o);
-      o=tmp_o;
-    }
-
-  /* Return */
-  return o;
-}
diff --git a/lib/arithmetic-bitand.c b/lib/arithmetic-bitand.c
new file mode 100644
index 0000000..d307285
--- /dev/null
+++ b/lib/arithmetic-bitand.c
@@ -0,0 +1,58 @@
+/*********************************************************************
+Arithmetic operations on data structures.
+This is part of GNU Astronomy Utilities (Gnuastro) package.
+
+Original author:
+     Mohammad Akhlaghi <address@hidden>
+Contributing author(s):
+Copyright (C) 2016, Free Software Foundation, Inc.
+
+Gnuastro is free software: you can redistribute it and/or modify it
+under the terms of the GNU General Public License as published by the
+Free Software Foundation, either version 3 of the License, or (at your
+option) any later version.
+
+Gnuastro is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
+**********************************************************************/
+#include <config.h>
+
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
+
+#include <gnuastro/blank.h>
+
+#include <gnuastro-internal/arithmetic-binary.h>
+
+
+
+
+
+
+/* `BINARY_SET_LT_INT' is defined in `arithmetic-binary.h'. As you see
+   there, this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_bitand(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  /* A small sanity check. */
+  if( l->type==GAL_TYPE_FLOAT32 || l->type==GAL_TYPE_FLOAT64
+      || r->type==GAL_TYPE_FLOAT32 || r->type==GAL_TYPE_FLOAT64 )
+      error(EXIT_FAILURE, 0, "%s: the bitand operator can only work on "
+            "integer type operands", __func__);
+
+  /* Do the operation. */
+  BINARY_SET_LT_INT( ( o->type==l->type
+                       ? ARITHMETIC_BINARY_OUT_TYPE_LEFT
+                       : ARITHMETIC_BINARY_OUT_TYPE_RIGHT ), & );
+}
diff --git a/lib/arithmetic-bitlsh.c b/lib/arithmetic-bitlsh.c
new file mode 100644
index 0000000..85538a1
--- /dev/null
+++ b/lib/arithmetic-bitlsh.c
@@ -0,0 +1,58 @@
+/*********************************************************************
+Arithmetic operations on data structures.
+This is part of GNU Astronomy Utilities (Gnuastro) package.
+
+Original author:
+     Mohammad Akhlaghi <address@hidden>
+Contributing author(s):
+Copyright (C) 2016, Free Software Foundation, Inc.
+
+Gnuastro is free software: you can redistribute it and/or modify it
+under the terms of the GNU General Public License as published by the
+Free Software Foundation, either version 3 of the License, or (at your
+option) any later version.
+
+Gnuastro is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
+**********************************************************************/
+#include <config.h>
+
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
+
+#include <gnuastro/blank.h>
+
+#include <gnuastro-internal/arithmetic-binary.h>
+
+
+
+
+
+
+/* `BINARY_SET_LT_INT' is defined in `arithmetic-binary.h'. As you see
+   there, this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_bitlsh(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  /* A small sanity check. */
+  if( l->type==GAL_TYPE_FLOAT32 || l->type==GAL_TYPE_FLOAT64
+      || r->type==GAL_TYPE_FLOAT32 || r->type==GAL_TYPE_FLOAT64 )
+      error(EXIT_FAILURE, 0, "%s: the bitlsh operator can only work on "
+            "integer type operands", __func__);
+
+  /* Do the operation. */
+  BINARY_SET_LT_INT( ( o->type==l->type
+                       ? ARITHMETIC_BINARY_OUT_TYPE_LEFT
+                       : ARITHMETIC_BINARY_OUT_TYPE_RIGHT ), << );
+}
diff --git a/lib/arithmetic-bitor.c b/lib/arithmetic-bitor.c
new file mode 100644
index 0000000..84ccea6
--- /dev/null
+++ b/lib/arithmetic-bitor.c
@@ -0,0 +1,58 @@
+/*********************************************************************
+Arithmetic operations on data structures.
+This is part of GNU Astronomy Utilities (Gnuastro) package.
+
+Original author:
+     Mohammad Akhlaghi <address@hidden>
+Contributing author(s):
+Copyright (C) 2016, Free Software Foundation, Inc.
+
+Gnuastro is free software: you can redistribute it and/or modify it
+under the terms of the GNU General Public License as published by the
+Free Software Foundation, either version 3 of the License, or (at your
+option) any later version.
+
+Gnuastro is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
+**********************************************************************/
+#include <config.h>
+
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
+
+#include <gnuastro/blank.h>
+
+#include <gnuastro-internal/arithmetic-binary.h>
+
+
+
+
+
+
+/* `BINARY_SET_LT_INT' is defined in `arithmetic-binary.h'. As you see
+   there, this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_bitor(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  /* A small sanity check. */
+  if( l->type==GAL_TYPE_FLOAT32 || l->type==GAL_TYPE_FLOAT64
+      || r->type==GAL_TYPE_FLOAT32 || r->type==GAL_TYPE_FLOAT64 )
+      error(EXIT_FAILURE, 0, "%s: the bitor operator can only work on "
+            "integer type operands", __func__);
+
+  /* Do the operation. */
+  BINARY_SET_LT_INT( ( o->type==l->type
+                       ? ARITHMETIC_BINARY_OUT_TYPE_LEFT
+                       : ARITHMETIC_BINARY_OUT_TYPE_RIGHT ), | );
+}
diff --git a/lib/arithmetic-bitrsh.c b/lib/arithmetic-bitrsh.c
new file mode 100644
index 0000000..a4f7e84
--- /dev/null
+++ b/lib/arithmetic-bitrsh.c
@@ -0,0 +1,58 @@
+/*********************************************************************
+Arithmetic operations on data structures.
+This is part of GNU Astronomy Utilities (Gnuastro) package.
+
+Original author:
+     Mohammad Akhlaghi <address@hidden>
+Contributing author(s):
+Copyright (C) 2016, Free Software Foundation, Inc.
+
+Gnuastro is free software: you can redistribute it and/or modify it
+under the terms of the GNU General Public License as published by the
+Free Software Foundation, either version 3 of the License, or (at your
+option) any later version.
+
+Gnuastro is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
+**********************************************************************/
+#include <config.h>
+
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
+
+#include <gnuastro/blank.h>
+
+#include <gnuastro-internal/arithmetic-binary.h>
+
+
+
+
+
+
+/* `BINARY_SET_LT_INT' is defined in `arithmetic-binary.h'. As you see
+   there, this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_bitrsh(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  /* A small sanity check. */
+  if( l->type==GAL_TYPE_FLOAT32 || l->type==GAL_TYPE_FLOAT64
+      || r->type==GAL_TYPE_FLOAT32 || r->type==GAL_TYPE_FLOAT64 )
+      error(EXIT_FAILURE, 0, "%s: the bitrsh operator can only work on "
+            "integer type operands", __func__);
+
+  /* Do the operation. */
+  BINARY_SET_LT_INT( ( o->type==l->type
+                       ? ARITHMETIC_BINARY_OUT_TYPE_LEFT
+                       : ARITHMETIC_BINARY_OUT_TYPE_RIGHT ), >> );
+}
diff --git a/lib/arithmetic-bitxor.c b/lib/arithmetic-bitxor.c
new file mode 100644
index 0000000..7df2f75
--- /dev/null
+++ b/lib/arithmetic-bitxor.c
@@ -0,0 +1,58 @@
+/*********************************************************************
+Arithmetic operations on data structures.
+This is part of GNU Astronomy Utilities (Gnuastro) package.
+
+Original author:
+     Mohammad Akhlaghi <address@hidden>
+Contributing author(s):
+Copyright (C) 2016, Free Software Foundation, Inc.
+
+Gnuastro is free software: you can redistribute it and/or modify it
+under the terms of the GNU General Public License as published by the
+Free Software Foundation, either version 3 of the License, or (at your
+option) any later version.
+
+Gnuastro is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
+**********************************************************************/
+#include <config.h>
+
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
+
+#include <gnuastro/blank.h>
+
+#include <gnuastro-internal/arithmetic-binary.h>
+
+
+
+
+
+
+/* `BINARY_SET_LT_INT' is defined in `arithmetic-binary.h'. As you see
+   there, this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_bitxor(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  /* A small sanity check. */
+  if( l->type==GAL_TYPE_FLOAT32 || l->type==GAL_TYPE_FLOAT64
+      || r->type==GAL_TYPE_FLOAT32 || r->type==GAL_TYPE_FLOAT64 )
+      error(EXIT_FAILURE, 0, "%s: the bitxor operator can only work on "
+            "integer type operands", __func__);
+
+  /* Do the operation. */
+  BINARY_SET_LT_INT( ( o->type==l->type
+                       ? ARITHMETIC_BINARY_OUT_TYPE_LEFT
+                       : ARITHMETIC_BINARY_OUT_TYPE_RIGHT ), ^ );
+}
diff --git a/lib/arithmetic-divide.c b/lib/arithmetic-divide.c
new file mode 100644
index 0000000..146bac3
--- /dev/null
+++ b/lib/arithmetic-divide.c
@@ -0,0 +1,53 @@
+/*********************************************************************
+Arithmetic operations on data structures.
+This is part of GNU Astronomy Utilities (Gnuastro) package.
+
+Original author:
+     Mohammad Akhlaghi <address@hidden>
+Contributing author(s):
+Copyright (C) 2016, Free Software Foundation, Inc.
+
+Gnuastro is free software: you can redistribute it and/or modify it
+under the terms of the GNU General Public License as published by the
+Free Software Foundation, either version 3 of the License, or (at your
+option) any later version.
+
+Gnuastro is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
+**********************************************************************/
+#include <config.h>
+
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
+
+#include <gnuastro/blank.h>
+
+#include <gnuastro-internal/arithmetic-binary.h>
+#include <gnuastro-internal/arithmetic-internal.h>
+
+
+
+
+
+/* `BINARY_SET_LT' is defined in `arithmetic-binary.h'. As you see there,
+   this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_divide(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  int checkblank=gal_arithmetic_binary_checkblank(l, r);
+
+  BINARY_SET_LT( ( o->type==l->type
+                   ? ARITHMETIC_BINARY_OUT_TYPE_LEFT
+                   : ARITHMETIC_BINARY_OUT_TYPE_RIGHT ), / );
+}
diff --git a/lib/gnuastro-internal/arithmetic-onlyint.h b/lib/arithmetic-eq.c
similarity index 50%
copy from lib/gnuastro-internal/arithmetic-onlyint.h
copy to lib/arithmetic-eq.c
index 2a8b3bd..8ce458a 100644
--- a/lib/gnuastro-internal/arithmetic-onlyint.h
+++ b/lib/arithmetic-eq.c
@@ -5,7 +5,7 @@ This is part of GNU Astronomy Utilities (Gnuastro) package.
 Original author:
      Mohammad Akhlaghi <address@hidden>
 Contributing author(s):
-Copyright (C) 2015, Free Software Foundation, Inc.
+Copyright (C) 2016, Free Software Foundation, Inc.
 
 Gnuastro is free software: you can redistribute it and/or modify it
 under the terms of the GNU General Public License as published by the
@@ -20,16 +20,32 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __GAL_ARITHMETIC_ONLYINT_H__
-#define __GAL_ARITHMETIC_ONLYINT_H__
+#include <config.h>
 
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
 
-gal_data_t *
-arithmetic_onlyint_binary(int operator, int flags, gal_data_t *lo,
-                          gal_data_t *ro);
+#include <gnuastro/blank.h>
 
-gal_data_t *
-arithmetic_onlyint_bitwise_not(int flags, gal_data_t *in);
+#include <gnuastro-internal/arithmetic-binary.h>
+#include <gnuastro-internal/arithmetic-internal.h>
 
 
-#endif
+
+
+
+/* `BINARY_SET_LT' is defined in `arithmetic-binary.h'. As you see there,
+   this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_eq(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  int checkblank=gal_arithmetic_binary_checkblank(l, r);
+
+  BINARY_SET_LT( ARITHMETIC_BINARY_OUT_TYPE_UINT8, == );
+}
diff --git a/lib/gnuastro-internal/arithmetic-onlyint.h b/lib/arithmetic-ge.c
similarity index 50%
copy from lib/gnuastro-internal/arithmetic-onlyint.h
copy to lib/arithmetic-ge.c
index 2a8b3bd..4f560af 100644
--- a/lib/gnuastro-internal/arithmetic-onlyint.h
+++ b/lib/arithmetic-ge.c
@@ -5,7 +5,7 @@ This is part of GNU Astronomy Utilities (Gnuastro) package.
 Original author:
      Mohammad Akhlaghi <address@hidden>
 Contributing author(s):
-Copyright (C) 2015, Free Software Foundation, Inc.
+Copyright (C) 2016, Free Software Foundation, Inc.
 
 Gnuastro is free software: you can redistribute it and/or modify it
 under the terms of the GNU General Public License as published by the
@@ -20,16 +20,32 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __GAL_ARITHMETIC_ONLYINT_H__
-#define __GAL_ARITHMETIC_ONLYINT_H__
+#include <config.h>
 
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
 
-gal_data_t *
-arithmetic_onlyint_binary(int operator, int flags, gal_data_t *lo,
-                          gal_data_t *ro);
+#include <gnuastro/blank.h>
 
-gal_data_t *
-arithmetic_onlyint_bitwise_not(int flags, gal_data_t *in);
+#include <gnuastro-internal/arithmetic-binary.h>
+#include <gnuastro-internal/arithmetic-internal.h>
 
 
-#endif
+
+
+
+/* `BINARY_SET_LT' is defined in `arithmetic-binary.h'. As you see there,
+   this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_ge(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  int checkblank=gal_arithmetic_binary_checkblank(l, r);
+
+  BINARY_SET_LT( ARITHMETIC_BINARY_OUT_TYPE_UINT8, >= );
+}
diff --git a/lib/gnuastro-internal/arithmetic-onlyint.h b/lib/arithmetic-gt.c
similarity index 50%
copy from lib/gnuastro-internal/arithmetic-onlyint.h
copy to lib/arithmetic-gt.c
index 2a8b3bd..f83deb7 100644
--- a/lib/gnuastro-internal/arithmetic-onlyint.h
+++ b/lib/arithmetic-gt.c
@@ -5,7 +5,7 @@ This is part of GNU Astronomy Utilities (Gnuastro) package.
 Original author:
      Mohammad Akhlaghi <address@hidden>
 Contributing author(s):
-Copyright (C) 2015, Free Software Foundation, Inc.
+Copyright (C) 2016, Free Software Foundation, Inc.
 
 Gnuastro is free software: you can redistribute it and/or modify it
 under the terms of the GNU General Public License as published by the
@@ -20,16 +20,32 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __GAL_ARITHMETIC_ONLYINT_H__
-#define __GAL_ARITHMETIC_ONLYINT_H__
+#include <config.h>
 
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
 
-gal_data_t *
-arithmetic_onlyint_binary(int operator, int flags, gal_data_t *lo,
-                          gal_data_t *ro);
+#include <gnuastro/blank.h>
 
-gal_data_t *
-arithmetic_onlyint_bitwise_not(int flags, gal_data_t *in);
+#include <gnuastro-internal/arithmetic-binary.h>
+#include <gnuastro-internal/arithmetic-internal.h>
 
 
-#endif
+
+
+
+/* `BINARY_SET_LT' is defined in `arithmetic-binary.h'. As you see there,
+   this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_gt(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  int checkblank=gal_arithmetic_binary_checkblank(l, r);
+
+  BINARY_SET_LT( ARITHMETIC_BINARY_OUT_TYPE_UINT8, > );
+}
diff --git a/lib/gnuastro-internal/arithmetic-onlyint.h b/lib/arithmetic-le.c
similarity index 50%
copy from lib/gnuastro-internal/arithmetic-onlyint.h
copy to lib/arithmetic-le.c
index 2a8b3bd..3d8f510 100644
--- a/lib/gnuastro-internal/arithmetic-onlyint.h
+++ b/lib/arithmetic-le.c
@@ -5,7 +5,7 @@ This is part of GNU Astronomy Utilities (Gnuastro) package.
 Original author:
      Mohammad Akhlaghi <address@hidden>
 Contributing author(s):
-Copyright (C) 2015, Free Software Foundation, Inc.
+Copyright (C) 2016, Free Software Foundation, Inc.
 
 Gnuastro is free software: you can redistribute it and/or modify it
 under the terms of the GNU General Public License as published by the
@@ -20,16 +20,32 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __GAL_ARITHMETIC_ONLYINT_H__
-#define __GAL_ARITHMETIC_ONLYINT_H__
+#include <config.h>
 
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
 
-gal_data_t *
-arithmetic_onlyint_binary(int operator, int flags, gal_data_t *lo,
-                          gal_data_t *ro);
+#include <gnuastro/blank.h>
 
-gal_data_t *
-arithmetic_onlyint_bitwise_not(int flags, gal_data_t *in);
+#include <gnuastro-internal/arithmetic-binary.h>
+#include <gnuastro-internal/arithmetic-internal.h>
 
 
-#endif
+
+
+
+/* `BINARY_SET_LT' is defined in `arithmetic-binary.h'. As you see there,
+   this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_le(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  int checkblank=gal_arithmetic_binary_checkblank(l, r);
+
+  BINARY_SET_LT( ARITHMETIC_BINARY_OUT_TYPE_UINT8, <= );
+}
diff --git a/lib/gnuastro-internal/arithmetic-onlyint.h b/lib/arithmetic-lt.c
similarity index 50%
copy from lib/gnuastro-internal/arithmetic-onlyint.h
copy to lib/arithmetic-lt.c
index 2a8b3bd..b99c874 100644
--- a/lib/gnuastro-internal/arithmetic-onlyint.h
+++ b/lib/arithmetic-lt.c
@@ -5,7 +5,7 @@ This is part of GNU Astronomy Utilities (Gnuastro) package.
 Original author:
      Mohammad Akhlaghi <address@hidden>
 Contributing author(s):
-Copyright (C) 2015, Free Software Foundation, Inc.
+Copyright (C) 2016, Free Software Foundation, Inc.
 
 Gnuastro is free software: you can redistribute it and/or modify it
 under the terms of the GNU General Public License as published by the
@@ -20,16 +20,32 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __GAL_ARITHMETIC_ONLYINT_H__
-#define __GAL_ARITHMETIC_ONLYINT_H__
+#include <config.h>
 
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
 
-gal_data_t *
-arithmetic_onlyint_binary(int operator, int flags, gal_data_t *lo,
-                          gal_data_t *ro);
+#include <gnuastro/blank.h>
 
-gal_data_t *
-arithmetic_onlyint_bitwise_not(int flags, gal_data_t *in);
+#include <gnuastro-internal/arithmetic-binary.h>
+#include <gnuastro-internal/arithmetic-internal.h>
 
 
-#endif
+
+
+
+/* `BINARY_SET_LT' is defined in `arithmetic-binary.h'. As you see there,
+   this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_lt(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  int checkblank=gal_arithmetic_binary_checkblank(l, r);
+
+  BINARY_SET_LT( ARITHMETIC_BINARY_OUT_TYPE_UINT8, < );
+}
diff --git a/lib/arithmetic-minus.c b/lib/arithmetic-minus.c
new file mode 100644
index 0000000..d88a0f6
--- /dev/null
+++ b/lib/arithmetic-minus.c
@@ -0,0 +1,53 @@
+/*********************************************************************
+Arithmetic operations on data structures.
+This is part of GNU Astronomy Utilities (Gnuastro) package.
+
+Original author:
+     Mohammad Akhlaghi <address@hidden>
+Contributing author(s):
+Copyright (C) 2016, Free Software Foundation, Inc.
+
+Gnuastro is free software: you can redistribute it and/or modify it
+under the terms of the GNU General Public License as published by the
+Free Software Foundation, either version 3 of the License, or (at your
+option) any later version.
+
+Gnuastro is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
+**********************************************************************/
+#include <config.h>
+
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
+
+#include <gnuastro/blank.h>
+
+#include <gnuastro-internal/arithmetic-binary.h>
+#include <gnuastro-internal/arithmetic-internal.h>
+
+
+
+
+
+/* `BINARY_SET_LT' is defined in `arithmetic-binary.h'. As you see there,
+   this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_minus(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  int checkblank=gal_arithmetic_binary_checkblank(l, r);
+
+  BINARY_SET_LT( ( o->type==l->type
+                   ? ARITHMETIC_BINARY_OUT_TYPE_LEFT
+                   : ARITHMETIC_BINARY_OUT_TYPE_RIGHT ), - );
+}
diff --git a/lib/arithmetic-modulo.c b/lib/arithmetic-modulo.c
new file mode 100644
index 0000000..f48db5e
--- /dev/null
+++ b/lib/arithmetic-modulo.c
@@ -0,0 +1,58 @@
+/*********************************************************************
+Arithmetic operations on data structures.
+This is part of GNU Astronomy Utilities (Gnuastro) package.
+
+Original author:
+     Mohammad Akhlaghi <address@hidden>
+Contributing author(s):
+Copyright (C) 2016, Free Software Foundation, Inc.
+
+Gnuastro is free software: you can redistribute it and/or modify it
+under the terms of the GNU General Public License as published by the
+Free Software Foundation, either version 3 of the License, or (at your
+option) any later version.
+
+Gnuastro is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
+**********************************************************************/
+#include <config.h>
+
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
+
+#include <gnuastro/blank.h>
+
+#include <gnuastro-internal/arithmetic-binary.h>
+
+
+
+
+
+
+/* `BINARY_SET_LT_INT' is defined in `arithmetic-binary.h'. As you see
+   there, this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_modulo(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  /* A small sanity check. */
+  if( l->type==GAL_TYPE_FLOAT32 || l->type==GAL_TYPE_FLOAT64
+      || r->type==GAL_TYPE_FLOAT32 || r->type==GAL_TYPE_FLOAT64 )
+      error(EXIT_FAILURE, 0, "%s: the modulo operator can only work on "
+            "integer type operands", __func__);
+
+  /* Do the operation. */
+  BINARY_SET_LT_INT( ( o->type==l->type
+                       ? ARITHMETIC_BINARY_OUT_TYPE_LEFT
+                       : ARITHMETIC_BINARY_OUT_TYPE_RIGHT ), % );
+}
diff --git a/lib/arithmetic-multiply.c b/lib/arithmetic-multiply.c
new file mode 100644
index 0000000..e060e36
--- /dev/null
+++ b/lib/arithmetic-multiply.c
@@ -0,0 +1,53 @@
+/*********************************************************************
+Arithmetic operations on data structures.
+This is part of GNU Astronomy Utilities (Gnuastro) package.
+
+Original author:
+     Mohammad Akhlaghi <address@hidden>
+Contributing author(s):
+Copyright (C) 2016, Free Software Foundation, Inc.
+
+Gnuastro is free software: you can redistribute it and/or modify it
+under the terms of the GNU General Public License as published by the
+Free Software Foundation, either version 3 of the License, or (at your
+option) any later version.
+
+Gnuastro is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
+**********************************************************************/
+#include <config.h>
+
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
+
+#include <gnuastro/blank.h>
+
+#include <gnuastro-internal/arithmetic-binary.h>
+#include <gnuastro-internal/arithmetic-internal.h>
+
+
+
+
+
+/* `BINARY_SET_LT' is defined in `arithmetic-binary.h'. As you see there,
+   this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_multiply(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  int checkblank=gal_arithmetic_binary_checkblank(l, r);
+
+  BINARY_SET_LT( ( o->type==l->type
+                   ? ARITHMETIC_BINARY_OUT_TYPE_LEFT
+                   : ARITHMETIC_BINARY_OUT_TYPE_RIGHT ), * );
+}
diff --git a/lib/gnuastro-internal/arithmetic-onlyint.h b/lib/arithmetic-ne.c
similarity index 50%
copy from lib/gnuastro-internal/arithmetic-onlyint.h
copy to lib/arithmetic-ne.c
index 2a8b3bd..3182b15 100644
--- a/lib/gnuastro-internal/arithmetic-onlyint.h
+++ b/lib/arithmetic-ne.c
@@ -5,7 +5,7 @@ This is part of GNU Astronomy Utilities (Gnuastro) package.
 Original author:
      Mohammad Akhlaghi <address@hidden>
 Contributing author(s):
-Copyright (C) 2015, Free Software Foundation, Inc.
+Copyright (C) 2016, Free Software Foundation, Inc.
 
 Gnuastro is free software: you can redistribute it and/or modify it
 under the terms of the GNU General Public License as published by the
@@ -20,16 +20,32 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __GAL_ARITHMETIC_ONLYINT_H__
-#define __GAL_ARITHMETIC_ONLYINT_H__
+#include <config.h>
 
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
 
-gal_data_t *
-arithmetic_onlyint_binary(int operator, int flags, gal_data_t *lo,
-                          gal_data_t *ro);
+#include <gnuastro/blank.h>
 
-gal_data_t *
-arithmetic_onlyint_bitwise_not(int flags, gal_data_t *in);
+#include <gnuastro-internal/arithmetic-binary.h>
+#include <gnuastro-internal/arithmetic-internal.h>
 
 
-#endif
+
+
+
+/* `BINARY_SET_LT' is defined in `arithmetic-binary.h'. As you see there,
+   this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_ne(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  int checkblank=gal_arithmetic_binary_checkblank(l, r);
+
+  BINARY_SET_LT( ARITHMETIC_BINARY_OUT_TYPE_UINT8, != );
+}
diff --git a/lib/arithmetic-onlyint.c b/lib/arithmetic-onlyint.c
deleted file mode 100644
index 9c60197..0000000
--- a/lib/arithmetic-onlyint.c
+++ /dev/null
@@ -1,492 +0,0 @@
-/*********************************************************************
-Arithmetic operations on data structures.
-This is part of GNU Astronomy Utilities (Gnuastro) package.
-
-Original author:
-     Mohammad Akhlaghi <address@hidden>
-Contributing author(s):
-Copyright (C) 2016, Free Software Foundation, Inc.
-
-Gnuastro is free software: you can redistribute it and/or modify it
-under the terms of the GNU General Public License as published by the
-Free Software Foundation, either version 3 of the License, or (at your
-option) any later version.
-
-Gnuastro is distributed in the hope that it will be useful, but
-WITHOUT ANY WARRANTY; without even the implied warranty of
-MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
-General Public License for more details.
-
-You should have received a copy of the GNU General Public License
-along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
-**********************************************************************/
-#include <config.h>
-
-#include <stdio.h>
-#include <errno.h>
-#include <error.h>
-#include <stdlib.h>
-
-#include <gnuastro/arithmetic.h>
-
-#include <gnuastro-internal/arithmetic-onlyint.h>
-#include <gnuastro-internal/arithmetic-internal.h>
-
-
-
-
-
-/************************************************************************/
-/*************            Native type macros            *****************/
-/************************************************************************/
-#if GAL_CONFIG_BIN_OP_UINT8 == 1
-#define BINOIN_LT_IS_UINT8                                         \
-  case GAL_TYPE_UINT8: BINOIN_LT_SET(uint8_t);          break;
-#define BINOIN_LT_SET_RT_IS_UINT8(LT)                              \
-  case GAL_TYPE_UINT8: BINOIN_RT_LT_SET(uint8_t, LT);   break;
-#else
-#define BINOIN_LT_IS_UINT8
-#define BINOIN_LT_SET_RT_IS_UINT8(LT)
-#endif
-
-
-
-
-
-#if GAL_CONFIG_BIN_OP_INT8 == 1
-#define BINOIN_LT_IS_INT8                                          \
-  case GAL_TYPE_INT8: BINOIN_LT_SET(int8_t);            break;
-#define BINOIN_LT_SET_RT_IS_INT8(LT)                               \
-  case GAL_TYPE_INT8: BINOIN_RT_LT_SET(int8_t, LT);     break;
-#else
-#define BINOIN_LT_IS_INT8
-#define BINOIN_LT_SET_RT_IS_INT8(LT)
-#endif
-
-
-
-
-
-#if GAL_CONFIG_BIN_OP_UINT16 == 1
-#define BINOIN_LT_IS_UINT16                                        \
-  case GAL_TYPE_UINT16: BINOIN_LT_SET(uint16_t);        break;
-#define BINOIN_LT_SET_RT_IS_UINT16(LT)                             \
-  case GAL_TYPE_UINT16: BINOIN_RT_LT_SET(uint16_t, LT); break;
-#else
-#define BINOIN_LT_IS_UINT16
-#define BINOIN_LT_SET_RT_IS_UINT16(LT)
-#endif
-
-
-
-
-
-#if GAL_CONFIG_BIN_OP_INT16 == 1
-#define BINOIN_LT_IS_INT16                                         \
-  case GAL_TYPE_INT16: BINOIN_LT_SET(int16_t);          break;
-#define BINOIN_LT_SET_RT_IS_INT16(LT)                              \
-  case GAL_TYPE_INT16: BINOIN_RT_LT_SET(int16_t, LT);   break;
-#else
-#define BINOIN_LT_IS_INT16
-#define BINOIN_LT_SET_RT_IS_INT16(LT)
-#endif
-
-
-
-
-
-#if GAL_CONFIG_BIN_OP_UINT32 == 1
-#define BINOIN_LT_IS_UINT32                                        \
-  case GAL_TYPE_UINT32: BINOIN_LT_SET(uint32_t);        break;
-#define BINOIN_LT_SET_RT_IS_UINT32(LT)                             \
-  case GAL_TYPE_UINT32: BINOIN_RT_LT_SET(uint32_t, LT); break;
-#else
-#define BINOIN_LT_IS_UINT32
-#define BINOIN_LT_SET_RT_IS_UINT32(LT)
-#endif
-
-
-
-
-
-#if GAL_CONFIG_BIN_OP_INT32 == 1
-#define BINOIN_LT_IS_INT32                                         \
-  case GAL_TYPE_INT32: BINOIN_LT_SET(int32_t);          break;
-#define BINOIN_LT_SET_RT_IS_INT32(LT)                              \
-  case GAL_TYPE_INT32: BINOIN_RT_LT_SET(int32_t, LT);   break;
-#else
-#define BINOIN_LT_IS_INT32
-#define BINOIN_LT_SET_RT_IS_INT32(LT)
-#endif
-
-
-
-
-
-#if GAL_CONFIG_BIN_OP_UINT64 == 1
-#define BINOIN_LT_IS_UINT64                                        \
-  case GAL_TYPE_UINT64: BINOIN_LT_SET(uint64_t);        break;
-#define BINOIN_LT_SET_RT_IS_UINT64(LT)                             \
-  case GAL_TYPE_UINT64: BINOIN_RT_LT_SET(uint64_t, LT); break;
-#else
-#define BINOIN_LT_IS_UINT64
-#define BINOIN_LT_SET_RT_IS_UINT64(LT)
-#endif
-
-
-
-
-
-#if GAL_CONFIG_BIN_OP_INT64 == 1
-#define BINOIN_LT_IS_INT64                                         \
-  case GAL_TYPE_INT64: BINOIN_LT_SET(int64_t);          break;
-#define BINOIN_LT_SET_RT_IS_INT64(LT)                              \
-  case GAL_TYPE_INT64: BINOIN_RT_LT_SET(int64_t, LT);   break;
-#else
-#define BINOIN_LT_IS_INT64
-#define BINOIN_LT_SET_RT_IS_INT64(LT)
-#endif
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-/************************************************************************/
-/*************              High level macros           *****************/
-/************************************************************************/
-/* Final step to be used by all operators and all types. */
-#define BINOIN_OP_OT_RT_LT_SET(OP, OT, RT, LT) {                   \
-    LT *la=l->array;                                               \
-    RT *ra=r->array;                                               \
-    OT *oa=o->array, *of=oa + o->size;                             \
-    if(l->size==r->size) do *oa = *la++ OP *ra++; while(++oa<of);  \
-    else if(l->size==1)  do *oa = *la   OP *ra++; while(++oa<of);  \
-    else                 do *oa = *la++ OP *ra;   while(++oa<of);  \
-  }
-
-
-
-
-
-/* For operators whose type may be any of the given inputs. */
-#define BINOIN_OP_RT_LT_SET(OP, RT, LT)                            \
-  if(o->type==l->type)                                             \
-    BINOIN_OP_OT_RT_LT_SET(OP, LT, RT, LT)                         \
-  else                                                             \
-    BINOIN_OP_OT_RT_LT_SET(OP, RT, RT, LT)
-
-
-
-
-
-/* Left and right types set, choose what to do based on operator. */
-#define BINOIN_RT_LT_SET(RT, LT)                                   \
-  switch(operator)                                                 \
-    {                                                              \
-    case GAL_ARITHMETIC_OP_MODULO:                                 \
-      BINOIN_OP_RT_LT_SET(%, RT, LT);                              \
-      break;                                                       \
-    case GAL_ARITHMETIC_OP_BITAND:                                 \
-      BINOIN_OP_RT_LT_SET(&, RT, LT);                              \
-      break;                                                       \
-    case GAL_ARITHMETIC_OP_BITOR:                                  \
-      BINOIN_OP_RT_LT_SET(|, RT, LT);                              \
-      break;                                                       \
-    case GAL_ARITHMETIC_OP_BITXOR:                                 \
-      BINOIN_OP_RT_LT_SET(^, RT, LT);                              \
-      break;                                                       \
-    case GAL_ARITHMETIC_OP_BITLSH:                                 \
-      BINOIN_OP_RT_LT_SET(<<, RT, LT);                             \
-      break;                                                       \
-    case GAL_ARITHMETIC_OP_BITRSH:                                 \
-      BINOIN_OP_RT_LT_SET(>>, RT, LT);                             \
-      break;                                                       \
-    default:                                                       \
-      error(EXIT_FAILURE, 0, "%s: operator code %d not recognized",\
-            "BINOIN_RT_LT_SET", operator);                         \
-    }
-
-
-
-
-
-
-/* Left operand type set, see what the right operand type is. */
-#define BINOIN_LT_SET(LT)                                          \
-  switch(r->type)                                                  \
-    {                                                              \
-      BINOIN_LT_SET_RT_IS_UINT8(LT);                               \
-      BINOIN_LT_SET_RT_IS_INT8(LT);                                \
-      BINOIN_LT_SET_RT_IS_UINT16(LT);                              \
-      BINOIN_LT_SET_RT_IS_INT16(LT);                               \
-      BINOIN_LT_SET_RT_IS_UINT32(LT);                              \
-      BINOIN_LT_SET_RT_IS_INT32(LT);                               \
-      BINOIN_LT_SET_RT_IS_UINT64(LT);                              \
-      BINOIN_LT_SET_RT_IS_INT64(LT);                               \
-    default:                                                       \
-      error(EXIT_FAILURE, 0, "%s: type code %d not recognized",    \
-            "BINOIN_LT_SET", r->type);                             \
-    }
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-/************************************************************************/
-/*************              Top level function          *****************/
-/************************************************************************/
-gal_data_t *
-arithmetic_onlyint_binary(int operator, int flags, gal_data_t *lo,
-                          gal_data_t *ro)
-{
-  /* Read the variable arguments. `lo' and `ro' keep the original data, in
-     case their type isn't built (based on configure options are configure
-     time). */
-  int otype, final_otype;
-  size_t out_size, minmapsize;
-  gal_data_t *l, *r, *o=NULL, *tmp_o;
-  char *opstring=gal_arithmetic_operator_string(operator);
-
-
-  /* Simple sanity check on the input sizes and types */
-  if( !( (flags & GAL_ARITHMETIC_NUMOK) && (lo->size==1 || ro->size==1))
-      && gal_data_dsize_is_different(lo, ro) )
-    error(EXIT_FAILURE, 0, "%s: the non-number inputs to %s don't have the "
-          "same dimension/size", __func__, opstring);
-
-  if( lo->type==GAL_TYPE_FLOAT32 || lo->type==GAL_TYPE_FLOAT64
-      || ro->type==GAL_TYPE_FLOAT32 || ro->type==GAL_TYPE_FLOAT64 )
-      error(EXIT_FAILURE, 0, "%s: the %s operator can only work on integer "
-            "type operands", __func__, opstring);
-
-
-  /* Set the final output type (independent of which types are
-     compiled). These needs to be done before the call to
-     `gal_arithmetic_convert_to_compiled_type', because that function can
-     free the space of the original data structures, thus we will loose the
-     original data structure information. */
-  final_otype=gal_arithmetic_binary_out_type(operator, lo, ro);
-
-
-  /* Make sure the input arrays have one of the compiled types. From this
-     point on, until the cleaning up section of this function, we won't be
-     using the `lo' and `ro' pointers. */
-  l=gal_arithmetic_convert_to_compiled_type(lo, flags);
-  r=gal_arithmetic_convert_to_compiled_type(ro, flags);
-
-
-  /* Sanity check: see if the compiled type is actually an integer. */
-  if( l->type>=GAL_TYPE_FLOAT32 || r->type>=GAL_TYPE_FLOAT32 )
-    error(EXIT_FAILURE, 0, "no larger integer compiled type. The `%s' "
-          "operator can only work on integer types. The left and right "
-          "operands had types `%s' and `%s'.\n\nYou can use the "
-          "`--enable-bin-op-XXXX' at configure time to compile a larger "
-          "type (note that unsigned types are considered to be larger than "
-          "signed ones). You can run the following command for more "
-          "information on these options (press the `SPACE' key to go down "
-          "and `q' to return to the command-line):\n\n"
-          "    $ info gnuastro \"Gnuastro configure options\"\n",
-          gal_arithmetic_operator_string(operator),
-          gal_type_name(lo->type, 1), gal_type_name(ro->type, 1));
-
-  /* Set the output type. */
-  otype=gal_type_out(l->type, r->type);
-
-
-  /* Set the output sizes. */
-  minmapsize = ( l->minmapsize < r->minmapsize
-                 ? l->minmapsize : r->minmapsize );
-  out_size = l->size > r->size ? l->size : r->size;
-
-
-  /* If we want inplace output, set the output pointer to one input. Note
-     that the output type can be different from both inputs.  */
-  if(flags & GAL_ARITHMETIC_INPLACE)
-    {
-      if     (l->type==otype && out_size==l->size)   o = l;
-      else if(r->type==otype && out_size==r->size)   o = r;
-    }
-
-
-  /* If the output pointer was not set for any reason, allocate it. For
-     `mmapsize', note that since its `size_t', it will always be
-     Positive. The `-1' that is recommended to give when you want the value
-     in RAM is actually the largest possible memory location. So we just
-     have to choose the smaller minmapsize of the two to decide if the
-     output array should be in RAM or not. */
-  if(o==NULL)
-    o = gal_data_alloc(NULL, otype,
-                       l->size>1 ? l->ndim  : r->ndim,
-                       l->size>1 ? l->dsize : r->dsize,
-                       l->size>1 ? l->wcs   : r->wcs,
-                       0, minmapsize, NULL, NULL, NULL );
-
-
-  /* Start setting the operator and operands. */
-  switch(l->type)
-    {
-      BINOIN_LT_IS_UINT8;
-      BINOIN_LT_IS_INT8;
-      BINOIN_LT_IS_UINT16;
-      BINOIN_LT_IS_INT16;
-      BINOIN_LT_IS_UINT32;
-      BINOIN_LT_IS_INT32;
-      BINOIN_LT_IS_UINT64;
-      BINOIN_LT_IS_INT64;
-    default:
-      error(EXIT_FAILURE, 0, "%s: type code %d not recognized",
-            __func__, l->type);
-    }
-
-
-  /* Clean up. Note that if the input arrays can be freed, and any of right
-     or left arrays needed conversion, `BINOIN_CONVERT_TO_COMPILED_TYPE'
-     has already freed the input arrays, so only `r' and `l' need
-     freeing. Alternatively, when the inputs shouldn't be freed, the only
-     allocated spaces are the `r' and `l' arrays if their types weren't
-     compiled for binary operations, we can tell this from the pointers: if
-     they are different from the original pointers, they were allocated. */
-  if(flags & GAL_ARITHMETIC_FREE)
-    {
-      if     (o==l)       gal_data_free(r);
-      else if(o==r)       gal_data_free(l);
-      else              { gal_data_free(l); gal_data_free(r); }
-    }
-  else
-    {
-      if(l!=lo)           gal_data_free(l);
-      if(r!=ro)           gal_data_free(r);
-    }
-
-
-  /* The type of the output dataset (`o->type') was chosen from `l' and `r'
-     (copies of the orignal operands but in a compiled type, not
-     necessarily the original `lo' and `ro' data structures). So we need to
-     to get the final type based on the original operands and check if the
-     final output needs changing.
-
-     IMPORTANT: This has to be done after (possibly) freeing the left and
-     right operands because this step can change the `o' pointer which
-     they may depend on (when working in-place). */
-  if( o->type != final_otype )
-    {
-      tmp_o=gal_data_copy_to_new_type(o, final_otype);
-      gal_data_free(o);
-      o=tmp_o;
-    }
-
-
-  /* Return */
-  return o;
-}
-
-
-
-
-
-gal_data_t *
-arithmetic_onlyint_bitwise_not(int flags, gal_data_t *in)
-{
-  gal_data_t *o;
-  uint8_t    *iu8  = in->array,  *iu8f  = iu8  + in->size,   *ou8;
-  int8_t     *ii8  = in->array,  *ii8f  = ii8  + in->size,   *oi8;
-  uint16_t   *iu16 = in->array,  *iu16f = iu16 + in->size,   *ou16;
-  int16_t    *ii16 = in->array,  *ii16f = ii16 + in->size,   *oi16;
-  uint32_t   *iu32 = in->array,  *iu32f = iu32 + in->size,   *ou32;
-  int32_t    *ii32 = in->array,  *ii32f = ii32 + in->size,   *oi32;
-  uint64_t   *iu64 = in->array,  *iu64f = iu64 + in->size,   *ou64;
-  int64_t    *ii64 = in->array,  *ii64f = ii64 + in->size,   *oi64;
-
-  /* Check the type */
-  switch(in->type)
-    {
-    case GAL_TYPE_FLOAT32:
-    case GAL_TYPE_FLOAT64:
-      error(EXIT_FAILURE, 0, "%s: bitwise not (one's complement) "
-            "operator can only work on integer types", __func__);
-    }
-
-  /* If we want inplace output, set the output pointer to the input
-     pointer, for every pixel, the operation will be independent. */
-  if(flags & GAL_ARITHMETIC_INPLACE)
-    o = in;
-  else
-    o = gal_data_alloc(NULL, in->type, in->ndim, in->dsize, in->wcs,
-                       0, in->minmapsize, NULL, NULL, NULL);
-
-  /* Start setting the types. */
-  switch(in->type)
-    {
-    case GAL_TYPE_UINT8:
-      ou8=o->array;   do  *ou8++ = ~(*iu8++);    while(iu8<iu8f);     break;
-
-    case GAL_TYPE_INT8:
-      oi8=o->array;   do  *oi8++ = ~(*ii8++);    while(ii8<ii8f);     break;
-
-    case GAL_TYPE_UINT16:
-      ou16=o->array;  do *ou16++ = ~(*iu16++);   while(iu16<iu16f);   break;
-
-    case GAL_TYPE_INT16:
-      oi16=o->array;  do *oi16++ = ~(*ii16++);   while(ii16<ii16f);   break;
-
-    case GAL_TYPE_UINT32:
-      ou32=o->array;  do *ou32++ = ~(*iu32++);   while(iu32<iu32f);   break;
-
-    case GAL_TYPE_INT32:
-      oi32=o->array;  do *oi32++ = ~(*ii32++);   while(ii32<ii32f);   break;
-
-    case GAL_TYPE_UINT64:
-      ou64=o->array;  do *ou64++ = ~(*iu64++);   while(iu64<iu64f);   break;
-
-    case GAL_TYPE_INT64:
-      oi64=o->array;  do *oi64++ = ~(*ii64++);   while(ii64<ii64f);   break;
-
-    default:
-      error(EXIT_FAILURE, 0, "%s: type code %d not recognized",
-            __func__, in->type);
-    }
-
-
-  /* Clean up. Note that if the input arrays can be freed, and any of right
-     or left arrays needed conversion, `UNIFUNC_CONVERT_TO_COMPILED_TYPE'
-     has already freed the input arrays, and we only have `r' and `l'
-     allocated in any case. Alternatively, when the inputs shouldn't be
-     freed, the only allocated spaces are the `r' and `l' arrays if their
-     types weren't compiled for binary operations, we can tell this from
-     the pointers: if they are different from the original pointers, they
-     were allocated. */
-  if( (flags & GAL_ARITHMETIC_FREE) && o!=in)
-    gal_data_free(in);
-
-  /* Return */
-  return o;
-}
diff --git a/lib/gnuastro-internal/arithmetic-onlyint.h b/lib/arithmetic-or.c
similarity index 50%
copy from lib/gnuastro-internal/arithmetic-onlyint.h
copy to lib/arithmetic-or.c
index 2a8b3bd..57be198 100644
--- a/lib/gnuastro-internal/arithmetic-onlyint.h
+++ b/lib/arithmetic-or.c
@@ -5,7 +5,7 @@ This is part of GNU Astronomy Utilities (Gnuastro) package.
 Original author:
      Mohammad Akhlaghi <address@hidden>
 Contributing author(s):
-Copyright (C) 2015, Free Software Foundation, Inc.
+Copyright (C) 2016, Free Software Foundation, Inc.
 
 Gnuastro is free software: you can redistribute it and/or modify it
 under the terms of the GNU General Public License as published by the
@@ -20,16 +20,32 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __GAL_ARITHMETIC_ONLYINT_H__
-#define __GAL_ARITHMETIC_ONLYINT_H__
+#include <config.h>
 
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
 
-gal_data_t *
-arithmetic_onlyint_binary(int operator, int flags, gal_data_t *lo,
-                          gal_data_t *ro);
+#include <gnuastro/blank.h>
 
-gal_data_t *
-arithmetic_onlyint_bitwise_not(int flags, gal_data_t *in);
+#include <gnuastro-internal/arithmetic-binary.h>
+#include <gnuastro-internal/arithmetic-internal.h>
 
 
-#endif
+
+
+
+/* `BINARY_SET_LT' is defined in `arithmetic-binary.h'. As you see there,
+   this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_or(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  int checkblank=gal_arithmetic_binary_checkblank(l, r);
+
+  BINARY_SET_LT( ARITHMETIC_BINARY_OUT_TYPE_UINT8, || );
+}
diff --git a/lib/arithmetic-plus.c b/lib/arithmetic-plus.c
new file mode 100644
index 0000000..22c7caa
--- /dev/null
+++ b/lib/arithmetic-plus.c
@@ -0,0 +1,53 @@
+/*********************************************************************
+Arithmetic operations on data structures.
+This is part of GNU Astronomy Utilities (Gnuastro) package.
+
+Original author:
+     Mohammad Akhlaghi <address@hidden>
+Contributing author(s):
+Copyright (C) 2016, Free Software Foundation, Inc.
+
+Gnuastro is free software: you can redistribute it and/or modify it
+under the terms of the GNU General Public License as published by the
+Free Software Foundation, either version 3 of the License, or (at your
+option) any later version.
+
+Gnuastro is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
+**********************************************************************/
+#include <config.h>
+
+#include <stdio.h>
+#include <errno.h>
+#include <error.h>
+#include <stdlib.h>
+
+#include <gnuastro/blank.h>
+
+#include <gnuastro-internal/arithmetic-binary.h>
+#include <gnuastro-internal/arithmetic-internal.h>
+
+
+
+
+
+/* `BINARY_SET_LT' is defined in `arithmetic-binary.h'. As you see there,
+   this is a deep macro (calls other macros) to deal with different
+   types. This allows efficiency in processing (after compilation), but
+   compilation will be very slow. Therefore, for each operator we have
+   defined a separate `.c' file so they are built separately and when built
+   in parallel can be much faster than having them all in a single file. */
+void
+arithmetic_plus(gal_data_t *l, gal_data_t *r, gal_data_t *o)
+{
+  int checkblank=gal_arithmetic_binary_checkblank(l, r);
+
+  BINARY_SET_LT( ( o->type==l->type
+                   ? ARITHMETIC_BINARY_OUT_TYPE_LEFT
+                   : ARITHMETIC_BINARY_OUT_TYPE_RIGHT ), + );
+}
diff --git a/lib/arithmetic.c b/lib/arithmetic.c
index 307284e..878173c 100644
--- a/lib/arithmetic.c
+++ b/lib/arithmetic.c
@@ -33,10 +33,84 @@ along with Gnuastro. If not, see 
<http://www.gnu.org/licenses/>.
 #include <gnuastro/statistics.h>
 #include <gnuastro/arithmetic.h>
 
-#include <gnuastro-internal/arithmetic-binary.h>
-#include <gnuastro-internal/arithmetic-onlyint.h>
 #include <gnuastro-internal/arithmetic-internal.h>
 
+/* Headers for each binary operator. Since they heavily involve macros,
+   their compilation can be very large if they are in a single function and
+   file. So there is a separate C source and header file for each of these
+   functions.*/
+#include <gnuastro-internal/arithmetic-lt.h>
+#include <gnuastro-internal/arithmetic-le.h>
+#include <gnuastro-internal/arithmetic-gt.h>
+#include <gnuastro-internal/arithmetic-ge.h>
+#include <gnuastro-internal/arithmetic-eq.h>
+#include <gnuastro-internal/arithmetic-ne.h>
+#include <gnuastro-internal/arithmetic-or.h>
+#include <gnuastro-internal/arithmetic-and.h>
+#include <gnuastro-internal/arithmetic-plus.h>
+#include <gnuastro-internal/arithmetic-minus.h>
+#include <gnuastro-internal/arithmetic-bitor.h>
+#include <gnuastro-internal/arithmetic-bitand.h>
+#include <gnuastro-internal/arithmetic-bitxor.h>
+#include <gnuastro-internal/arithmetic-bitlsh.h>
+#include <gnuastro-internal/arithmetic-bitrsh.h>
+#include <gnuastro-internal/arithmetic-modulo.h>
+#include <gnuastro-internal/arithmetic-divide.h>
+#include <gnuastro-internal/arithmetic-multiply.h>
+
+
+
+
+
+
+
+
+
+
+/***********************************************************************/
+/***************             Internal checks              **************/
+/***********************************************************************/
+/* Some functions are only for a floating point operand, so if the input
+   isn't floating point, inform the user to change the type explicitly,
+   doing it implicitly/internally puts too much responsability on the
+   program. */
+static void
+arithmetic_check_float_input(gal_data_t *in, int operator, char *numstr)
+{
+  switch(in->type)
+    {
+    case GAL_TYPE_FLOAT32:
+    case GAL_TYPE_FLOAT64:
+      break;
+    default:
+      error(EXIT_FAILURE, 0, "the %s operator can only accept single or "
+            "double precision floating point numbers as its operand. The "
+            "%s operand has type %s. You can use the `float' or `double' "
+            "operators before this operator to explicitly convert to the "
+            "desired precision floating point type. If the operand was "
+            "originally a typed number (string of characters), add an `f' "
+            "after it so it is directly read into the proper precision "
+            "floating point number (based on the number of non-zero "
+            "decimals it has)", gal_arithmetic_operator_string(operator),
+            numstr, gal_type_name(in->type, 1));
+    }
+}
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
 
 
 
@@ -139,6 +213,81 @@ arithmetic_not(gal_data_t *data, int flags)
 
 
 
+/* Bitwise not operator. */
+static gal_data_t *
+arithmetic_bitwise_not(int flags, gal_data_t *in)
+{
+  gal_data_t *o;
+  uint8_t    *iu8  = in->array,  *iu8f  = iu8  + in->size,   *ou8;
+  int8_t     *ii8  = in->array,  *ii8f  = ii8  + in->size,   *oi8;
+  uint16_t   *iu16 = in->array,  *iu16f = iu16 + in->size,   *ou16;
+  int16_t    *ii16 = in->array,  *ii16f = ii16 + in->size,   *oi16;
+  uint32_t   *iu32 = in->array,  *iu32f = iu32 + in->size,   *ou32;
+  int32_t    *ii32 = in->array,  *ii32f = ii32 + in->size,   *oi32;
+  uint64_t   *iu64 = in->array,  *iu64f = iu64 + in->size,   *ou64;
+  int64_t    *ii64 = in->array,  *ii64f = ii64 + in->size,   *oi64;
+
+  /* Check the type */
+  switch(in->type)
+    {
+    case GAL_TYPE_FLOAT32:
+    case GAL_TYPE_FLOAT64:
+      error(EXIT_FAILURE, 0, "%s: bitwise not (one's complement) "
+            "operator can only work on integer types", __func__);
+    }
+
+  /* If we want inplace output, set the output pointer to the input
+     pointer, for every pixel, the operation will be independent. */
+  if(flags & GAL_ARITHMETIC_INPLACE)
+    o = in;
+  else
+    o = gal_data_alloc(NULL, in->type, in->ndim, in->dsize, in->wcs,
+                       0, in->minmapsize, NULL, NULL, NULL);
+
+  /* Start setting the types. */
+  switch(in->type)
+    {
+    case GAL_TYPE_UINT8:
+      ou8=o->array;   do  *ou8++ = ~(*iu8++);    while(iu8<iu8f);     break;
+
+    case GAL_TYPE_INT8:
+      oi8=o->array;   do  *oi8++ = ~(*ii8++);    while(ii8<ii8f);     break;
+
+    case GAL_TYPE_UINT16:
+      ou16=o->array;  do *ou16++ = ~(*iu16++);   while(iu16<iu16f);   break;
+
+    case GAL_TYPE_INT16:
+      oi16=o->array;  do *oi16++ = ~(*ii16++);   while(ii16<ii16f);   break;
+
+    case GAL_TYPE_UINT32:
+      ou32=o->array;  do *ou32++ = ~(*iu32++);   while(iu32<iu32f);   break;
+
+    case GAL_TYPE_INT32:
+      oi32=o->array;  do *oi32++ = ~(*ii32++);   while(ii32<ii32f);   break;
+
+    case GAL_TYPE_UINT64:
+      ou64=o->array;  do *ou64++ = ~(*iu64++);   while(iu64<iu64f);   break;
+
+    case GAL_TYPE_INT64:
+      oi64=o->array;  do *oi64++ = ~(*ii64++);   while(ii64<ii64f);   break;
+
+    default:
+      error(EXIT_FAILURE, 0, "%s: type code %d not recognized",
+            __func__, in->type);
+    }
+
+
+  /* Clean up (if necessary). */
+  if( (flags & GAL_ARITHMETIC_FREE) && o!=in)
+    gal_data_free(in);
+
+  /* Return */
+  return o;
+}
+
+
+
+
 
 /* We don't want to use the standard function for unary functions in the
    case of the absolute operator. This is because there are multiple
@@ -198,82 +347,11 @@ arithmetic_abs(int flags, gal_data_t *in)
 
 
 
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-/***********************************************************************/
-/***************          Checking functions              **************/
-/***********************************************************************/
-/* Some functions are only for a floating point operand, so if the input
-   isn't floating point, inform the user to change the type explicitly,
-   doing it implicitly/internally puts too much responsability on the
-   program. */
-static void
-arithmetic_check_float_input(gal_data_t *in, int operator, char *numstr)
-{
-  switch(in->type)
-    {
-    case GAL_TYPE_FLOAT32:
-    case GAL_TYPE_FLOAT64:
-      break;
-    default:
-      error(EXIT_FAILURE, 0, "the %s operator can only accept single or "
-            "double precision floating point numbers as its operand. The "
-            "%s operand has type %s. You can use the `float' or `double' "
-            "operators before this operator to explicitly convert to the "
-            "desired precision floating point type. If the operand was "
-            "originally a typed number (string of characters), add an `f' "
-            "after it so it is directly read into the proper precision "
-            "floating point number (based on the number of non-zero "
-            "decimals it has)", gal_arithmetic_operator_string(operator),
-            numstr, gal_type_name(in->type, 1));
-    }
-}
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-/***********************************************************************/
-/***************             Unary functions              **************/
-/***********************************************************************/
-
 #define UNIFUNC_RUN_FUNCTION_ON_ELEMENT(IT, OP){                        \
     IT *ia=in->array, *oa=o->array, *iaf=ia + in->size;                 \
     do *oa++ = OP(*ia++); while(ia<iaf);                                \
   }
 
-
-
-
-
 #define UNIARY_FUNCTION_ON_ELEMENT(OP)                                  \
   switch(in->type)                                                      \
     {                                                                   \
@@ -312,10 +390,6 @@ arithmetic_check_float_input(gal_data_t *in, int operator, 
char *numstr)
             "UNIARY_FUNCTION_ON_ELEMENT", in->type);                    \
     }
 
-
-
-
-
 static gal_data_t *
 arithmetic_unary_function(int operator, int flags, gal_data_t *in)
 {
@@ -369,7 +443,32 @@ arithmetic_unary_function(int operator, int flags, 
gal_data_t *in)
 
 
 
+/* Call functions in the `gnuastro/statistics' library. */
+static gal_data_t *
+arithmetic_from_statistics(int operator, int flags, gal_data_t *input)
+{
+  gal_data_t *out=NULL;
+  int ip=(flags & GAL_ARITHMETIC_INPLACE) || (flags & GAL_ARITHMETIC_FREE);
+
+  switch(operator)
+    {
+    case GAL_ARITHMETIC_OP_MINVAL:  out=gal_statistics_minimum(input); break;
+    case GAL_ARITHMETIC_OP_MAXVAL:  out=gal_statistics_maximum(input); break;
+    case GAL_ARITHMETIC_OP_NUMVAL:  out=gal_statistics_number(input);  break;
+    case GAL_ARITHMETIC_OP_SUMVAL:  out=gal_statistics_sum(input);     break;
+    case GAL_ARITHMETIC_OP_MEANVAL: out=gal_statistics_mean(input);    break;
+    case GAL_ARITHMETIC_OP_STDVAL:  out=gal_statistics_std(input);     break;
+    case GAL_ARITHMETIC_OP_MEDIANVAL:
+      out=gal_statistics_median(input, ip); break;
+    default:
+      error(EXIT_FAILURE, 0, "%s: operator code %d not recognized",
+            __func__, operator);
+    }
 
+  /* If the input is to be freed, then do so and return the output. */
+  if( flags & GAL_ARITHMETIC_FREE ) gal_data_free(input);
+  return out;
+}
 
 
 
@@ -384,207 +483,42 @@ arithmetic_unary_function(int operator, int flags, 
gal_data_t *in)
 
 
 
-/***********************************************************************/
-/***************            Binary functions              **************/
-/***********************************************************************/
 
 
-#define BINFUNC_RUN_FUNCTION(OT, RT, LT, OP){                           \
-    LT *la=l->array;                                                    \
-    RT *ra=r->array;                                                    \
-    OT *oa=o->array, *of=oa + o->size;                                  \
-    if(l->size==r->size) do *oa = OP(*la++, *ra++); while(++oa<of);     \
-    else if(l->size==1)  do *oa = OP(*la,   *ra++); while(++oa<of);     \
-    else                 do *oa = OP(*la++, *ra  ); while(++oa<of);     \
-  }
 
 
 
 
+/***********************************************************************/
+/***************                  Where                   **************/
+/***********************************************************************/
 
-#define BINFUNC_F_OPERATOR_LEFT_RIGHT_SET(RT, LT, OP)                   \
-  switch(o->type)                                                       \
-    {                                                                   \
-    case GAL_TYPE_FLOAT32:                                              \
-      BINFUNC_RUN_FUNCTION(float, RT, LT, OP);                          \
-      break;                                                            \
-    case GAL_TYPE_FLOAT64:                                              \
-      BINFUNC_RUN_FUNCTION(double, RT, LT, OP);                         \
-      break;                                                            \
-    default:                                                            \
-      error(EXIT_FAILURE, 0, "%s: type %d not recognized for o->type ", \
-            "BINFUNC_F_OPERATOR_LEFT_RIGHT_SET", o->type);              \
-    }
+/* When the `iftrue' dataset only has one element and the element is blank,
+   then it will be replaced with the blank value of the type of the output
+   data. */
+#define DO_WHERE_OPERATION(ITT, OT) {                                \
+    ITT *it=iftrue->array;                                           \
+    OT b, *o=out->array, *of=o+out->size;                            \
+    if(iftrue->size==1)                                              \
+      {                                                              \
+        if( gal_blank_present(iftrue, 0) )                           \
+          {                                                          \
+            gal_blank_write(&b, out->type);                          \
+            do { *o = *c++ ? b : *o;        } while(++o<of);         \
+          }                                                          \
+        else                                                         \
+          do   { *o = *c++ ? *it : *o;       } while(++o<of);        \
+      }                                                              \
+    else                                                             \
+      do       { *o = *c++ ? *it : *o; ++it; } while(++o<of);        \
+}
 
 
 
 
 
-#define BINFUNC_F_OPERATOR_LEFT_SET(LT, OP)                             \
-  switch(r->type)                                                       \
-    {                                                                   \
-    case GAL_TYPE_FLOAT32:                                              \
-      BINFUNC_F_OPERATOR_LEFT_RIGHT_SET(float, LT, OP);                 \
-      break;                                                            \
-    case GAL_TYPE_FLOAT64:                                              \
-      BINFUNC_F_OPERATOR_LEFT_RIGHT_SET(double, LT, OP);                \
-      break;                                                            \
-    default:                                                            \
-      error(EXIT_FAILURE, 0, "%s: type %d not recognized for r->type",  \
-            "BINFUNC_F_OPERATOR_LEFT_SET", r->type);                    \
-    }
-
-
-
-
-
-#define BINFUNC_F_OPERATOR_SET(OP)                                      \
-  switch(l->type)                                                       \
-    {                                                                   \
-    case GAL_TYPE_FLOAT32:                                              \
-      BINFUNC_F_OPERATOR_LEFT_SET(float, OP);                           \
-      break;                                                            \
-    case GAL_TYPE_FLOAT64:                                              \
-      BINFUNC_F_OPERATOR_LEFT_SET(double, OP);                          \
-      break;                                                            \
-    default:                                                            \
-      error(EXIT_FAILURE, 0, "%s: type %d not recognized for l->type",  \
-            "BINFUNC_F_OPERATOR_SET", l->type);                         \
-    }
-
-
-
-
-
-static gal_data_t *
-arithmetic_binary_function_flt(int operator, int flags, gal_data_t *l,
-                               gal_data_t *r)
-{
-  int final_otype;
-  gal_data_t *o=NULL;
-  size_t out_size, minmapsize;
-
-
-  /* Simple sanity check on the input sizes */
-  if( !( (flags & GAL_ARITHMETIC_NUMOK) && (l->size==1 || r->size==1))
-      && gal_data_dsize_is_different(l, r) )
-    error(EXIT_FAILURE, 0, "%s: the input datasets don't have the same "
-          "dimension/size", __func__);
-
-  /* Check for the types of the left and right operands. */
-  arithmetic_check_float_input(l, operator, "first");
-  arithmetic_check_float_input(r, operator, "second");
-
-  /* Set the output type. */
-  final_otype = gal_type_out(l->type, r->type);
-
-  /* Set the output sizes. */
-  minmapsize = ( l->minmapsize < r->minmapsize
-                 ? l->minmapsize : r->minmapsize );
-  out_size = l->size > r->size ? l->size : r->size;
-
-
-  /* If we want inplace output, set the output pointer to one input. Note
-     that the output type can be different from both inputs.  */
-  if(flags & GAL_ARITHMETIC_INPLACE)
-    {
-      if     (l->type==final_otype && out_size==l->size)   o = l;
-      else if(r->type==final_otype && out_size==r->size)   o = r;
-    }
-
-
-  /* If the output pointer was not set for any reason, allocate it. For
-     `mmapsize', note that since its `size_t', it will always be
-     Positive. The `-1' that is recommended to give when you want the value
-     in RAM is actually the largest possible memory location. So we just
-     have to choose the smaller minmapsize of the two to decide if the
-     output array should be in RAM or not. */
-  if(o==NULL)
-    o = gal_data_alloc(NULL, final_otype,
-                       l->size>1 ? l->ndim  : r->ndim,
-                       l->size>1 ? l->dsize : r->dsize,
-                       l->size>1 ? l->wcs : r->wcs, 0, minmapsize,
-                       NULL, NULL, NULL);
-
-
-  /* Start setting the operator and operands. */
-  switch(operator)
-    {
-    case GAL_ARITHMETIC_OP_POW:  BINFUNC_F_OPERATOR_SET( pow  ); break;
-    default:
-      error(EXIT_FAILURE, 0, "%s: operator code %d not recognized",
-            __func__, operator);
-    }
-
-
-  /* Clean up. Note that if the input arrays can be freed, and any of right
-     or left arrays needed conversion, `BINFUNC_CONVERT_TO_COMPILED_TYPE'
-     has already freed the input arrays, and we only have `r' and `l'
-     allocated in any case. Alternatively, when the inputs shouldn't be
-     freed, the only allocated spaces are the `r' and `l' arrays if their
-     types weren't compiled for binary operations, we can tell this from
-     the pointers: if they are different from the original pointers, they
-     were allocated. */
-  if(flags & GAL_ARITHMETIC_FREE)
-    {
-      if     (o==l)       gal_data_free(r);
-      else if(o==r)       gal_data_free(l);
-      else              { gal_data_free(l); gal_data_free(r); }
-    }
-
-  /* Return */
-  return o;
-}
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-/***********************************************************************/
-/***************                  Where                   **************/
-/***********************************************************************/
-
-/* When the `iftrue' dataset only has one element and the element is blank,
-   then it will be replaced with the blank value of the type of the output
-   data. */
-#define DO_WHERE_OPERATION(ITT, OT) {                                \
-    ITT *it=iftrue->array;                                           \
-    OT b, *o=out->array, *of=o+out->size;                            \
-    if(iftrue->size==1)                                              \
-      {                                                              \
-        if( gal_blank_present(iftrue, 0) )                           \
-          {                                                          \
-            gal_blank_write(&b, out->type);                          \
-            do { *o = *c++ ? b : *o;        } while(++o<of);         \
-          }                                                          \
-        else                                                         \
-          do   { *o = *c++ ? *it : *o;       } while(++o<of);        \
-      }                                                              \
-    else                                                             \
-      do       { *o = *c++ ? *it : *o; ++it; } while(++o<of);        \
-}
-
-
-
-
-
-#define WHERE_OUT_SET(OT)                                               \
-  switch(iftrue->type)                                                  \
+#define WHERE_OUT_SET(OT)                                               \
+  switch(iftrue->type)                                                  \
     {                                                                   \
     case GAL_TYPE_UINT8:    DO_WHERE_OPERATION( uint8_t,  OT);  break;  \
     case GAL_TYPE_INT8:     DO_WHERE_OPERATION( int8_t,   OT);  break;  \
@@ -1083,10 +1017,28 @@ arithmetic_multioperand(int operator, int flags, 
gal_data_t *list)
 
 
 /**********************************************************************/
-/****************      Compiled binary op types       *****************/
+/****************           Binary operators          *****************/
 /**********************************************************************/
+/* See if we should check for blanks. When both types are floats, blanks
+   don't need to be checked (the floating point standard will do the job
+   for us). It is also not necessary to check blanks in bitwise operators,
+   but bitwise operators have their own macro
+   (`BINARY_OP_INCR_OT_RT_LT_SET') which doesn' use `checkblanks'.*/
 int
-gal_arithmetic_binary_out_type(int operator, gal_data_t *l, gal_data_t *r)
+gal_arithmetic_binary_checkblank(gal_data_t *l, gal_data_t *r)
+{
+  return ((((l->type!=GAL_TYPE_FLOAT32    && l->type!=GAL_TYPE_FLOAT64)
+            || (r->type!=GAL_TYPE_FLOAT32 && r->type!=GAL_TYPE_FLOAT64))
+           && (gal_blank_present(l, 1) || gal_blank_present(r, 1)))
+          ? 1 : 0 );
+}
+
+
+
+
+
+static int
+arithmetic_binary_out_type(int operator, gal_data_t *l, gal_data_t *r)
 {
   switch(operator)
     {
@@ -1106,118 +1058,245 @@ gal_arithmetic_binary_out_type(int operator, 
gal_data_t *l, gal_data_t *r)
 
 
 
-static int
-arithmetic_nearest_compiled_type(int intype)
+static gal_data_t *
+arithmetic_binary(int operator, int flags, gal_data_t *l, gal_data_t *r)
 {
-  switch(intype)
+  /* Read the variable arguments. `lo' and `ro' keep the original data, in
+     case their type isn't built (based on configure options are configure
+     time). */
+  int32_t otype;
+  gal_data_t *o=NULL;
+  size_t out_size, minmapsize;
+
+
+  /* Simple sanity check on the input sizes */
+  if( !( (flags & GAL_ARITHMETIC_NUMOK) && (l->size==1 || r->size==1))
+      && gal_data_dsize_is_different(l, r) )
+    error(EXIT_FAILURE, 0, "%s: the non-number inputs to %s don't have the "
+          "same dimension/size", __func__,
+          gal_arithmetic_operator_string(operator));
+
+
+  /* Set the output type. For the comparison operators, the output type is
+     either 0 or 1, so we will set the output type to `unsigned char' for
+     efficient memory and CPU usage. Since the number of operators without
+     a fixed output type (like the conditionals) is less, by `default' we
+     will set the output type to `unsigned char', and if any of the other
+     operatrs are given, it will be chosen based on the input types.*/
+  otype=arithmetic_binary_out_type(operator, l, r);
+
+
+  /* Set the output sizes. */
+  minmapsize = ( l->minmapsize < r->minmapsize
+                 ? l->minmapsize : r->minmapsize );
+  out_size = l->size > r->size ? l->size : r->size;
+
+
+  /* If we want inplace output, set the output pointer to one input. Note
+     that the output type can be different from both inputs.  */
+  if(flags & GAL_ARITHMETIC_INPLACE)
     {
-    case GAL_TYPE_UINT8:
-      if(GAL_CONFIG_BIN_OP_UINT8) return GAL_TYPE_UINT8;
-      else
-        {
-          if     (GAL_CONFIG_BIN_OP_UINT16)   return GAL_TYPE_UINT16;
-          else if(GAL_CONFIG_BIN_OP_INT16)    return GAL_TYPE_INT16;
-          else if(GAL_CONFIG_BIN_OP_UINT32)   return GAL_TYPE_UINT32;
-          else if(GAL_CONFIG_BIN_OP_INT32)    return GAL_TYPE_INT32;
-          else if(GAL_CONFIG_BIN_OP_UINT64)   return GAL_TYPE_UINT64;
-          else if(GAL_CONFIG_BIN_OP_INT64)    return GAL_TYPE_INT64;
-          else if(GAL_CONFIG_BIN_OP_FLOAT32)  return GAL_TYPE_FLOAT32;
-          else if(GAL_CONFIG_BIN_OP_FLOAT64)  return GAL_TYPE_FLOAT64;
-        }
-      break;
+      if     (l->type==otype && out_size==l->size)   o = l;
+      else if(r->type==otype && out_size==r->size)   o = r;
+    }
 
-    case GAL_TYPE_INT8:
-      if(GAL_CONFIG_BIN_OP_INT8) return GAL_TYPE_INT8;
-      else
-        {
-          if     (GAL_CONFIG_BIN_OP_INT16)    return GAL_TYPE_INT16;
-          else if(GAL_CONFIG_BIN_OP_INT32)    return GAL_TYPE_INT32;
-          else if(GAL_CONFIG_BIN_OP_INT64)    return GAL_TYPE_INT64;
-          else if(GAL_CONFIG_BIN_OP_FLOAT32)  return GAL_TYPE_FLOAT32;
-          else if(GAL_CONFIG_BIN_OP_FLOAT64)  return GAL_TYPE_FLOAT64;
-        }
-      break;
 
-    case GAL_TYPE_UINT16:
-      if(GAL_CONFIG_BIN_OP_UINT16) return GAL_TYPE_UINT16;
-      else
-        {
-          if     (GAL_CONFIG_BIN_OP_UINT32)   return GAL_TYPE_UINT32;
-          else if(GAL_CONFIG_BIN_OP_INT32)    return GAL_TYPE_INT32;
-          else if(GAL_CONFIG_BIN_OP_UINT64)   return GAL_TYPE_UINT64;
-          else if(GAL_CONFIG_BIN_OP_INT64)    return GAL_TYPE_INT64;
-          else if(GAL_CONFIG_BIN_OP_FLOAT32)  return GAL_TYPE_FLOAT32;
-          else if(GAL_CONFIG_BIN_OP_FLOAT64)  return GAL_TYPE_FLOAT64;
-        }
-      break;
+  /* If the output pointer was not set above for any of the possible
+     reasons, allocate it. For `mmapsize', note that since its `size_t', it
+     will always be positive. The `-1' that is recommended to give when you
+     want the value in RAM is actually the largest possible memory
+     location. So we just have to choose the smaller minmapsize of the two
+     to decide if the output array should be in RAM or not. */
+  if(o==NULL)
+    o = gal_data_alloc(NULL, otype,
+                       l->size>1 ? l->ndim  : r->ndim,
+                       l->size>1 ? l->dsize : r->dsize,
+                       l->size>1 ? l->wcs   : r->wcs,
+                       0, minmapsize, NULL, NULL, NULL );
 
-    case GAL_TYPE_INT16:
-      if(GAL_CONFIG_BIN_OP_INT16) return GAL_TYPE_INT16;
-      else
-        {
-          if     (GAL_CONFIG_BIN_OP_INT32)    return GAL_TYPE_INT32;
-          else if(GAL_CONFIG_BIN_OP_INT64)    return GAL_TYPE_INT64;
-          else if(GAL_CONFIG_BIN_OP_FLOAT32)  return GAL_TYPE_FLOAT32;
-          else if(GAL_CONFIG_BIN_OP_FLOAT64)  return GAL_TYPE_FLOAT64;
-        }
-      break;
 
-    case GAL_TYPE_UINT32:
-      if(GAL_CONFIG_BIN_OP_UINT32) return GAL_TYPE_UINT32;
-      else
-        {
-          if     (GAL_CONFIG_BIN_OP_UINT64)   return GAL_TYPE_UINT64;
-          else if(GAL_CONFIG_BIN_OP_INT64)    return GAL_TYPE_INT64;
-          else if(GAL_CONFIG_BIN_OP_FLOAT32)  return GAL_TYPE_FLOAT32;
-          else if(GAL_CONFIG_BIN_OP_FLOAT64)  return GAL_TYPE_FLOAT64;
-        }
-      break;
+  /* Call the proper function for the operator. Since they heavily involve
+     macros, their compilation can be very large if they are in a single
+     function and file. So there is a separate C source and header file for
+     each of these functions. */
+  switch(operator)
+    {
+    case GAL_ARITHMETIC_OP_PLUS:     arithmetic_plus(l, r, o);     break;
+    case GAL_ARITHMETIC_OP_MINUS:    arithmetic_minus(l, r, o);    break;
+    case GAL_ARITHMETIC_OP_MULTIPLY: arithmetic_multiply(l, r, o); break;
+    case GAL_ARITHMETIC_OP_DIVIDE:   arithmetic_divide(l, r, o);   break;
+    case GAL_ARITHMETIC_OP_LT:       arithmetic_lt(l, r, o);       break;
+    case GAL_ARITHMETIC_OP_LE:       arithmetic_le(l, r, o);       break;
+    case GAL_ARITHMETIC_OP_GT:       arithmetic_gt(l, r, o);       break;
+    case GAL_ARITHMETIC_OP_GE:       arithmetic_ge(l, r, o);       break;
+    case GAL_ARITHMETIC_OP_EQ:       arithmetic_eq(l, r, o);       break;
+    case GAL_ARITHMETIC_OP_NE:       arithmetic_ne(l, r, o);       break;
+    case GAL_ARITHMETIC_OP_AND:      arithmetic_and(l, r, o);      break;
+    case GAL_ARITHMETIC_OP_OR:       arithmetic_or(l, r, o);       break;
+    case GAL_ARITHMETIC_OP_BITAND:   arithmetic_bitand(l, r, o);   break;
+    case GAL_ARITHMETIC_OP_BITOR:    arithmetic_bitor(l, r, o);    break;
+    case GAL_ARITHMETIC_OP_BITXOR:   arithmetic_bitxor(l, r, o);   break;
+    case GAL_ARITHMETIC_OP_BITLSH:   arithmetic_bitlsh(l, r, o);   break;
+    case GAL_ARITHMETIC_OP_BITRSH:   arithmetic_bitrsh(l, r, o);   break;
+    case GAL_ARITHMETIC_OP_MODULO:   arithmetic_modulo(l, r, o);   break;
+    default:
+      error(EXIT_FAILURE, 0, "%s: a bug! please contact us at %s to address "
+            "the problem. %d is not a valid operator code", __func__,
+            PACKAGE_BUGREPORT, operator);
+    }
 
-    case GAL_TYPE_INT32:
-      if(GAL_CONFIG_BIN_OP_INT32) return GAL_TYPE_INT32;
-      else
-        {
-          if     (GAL_CONFIG_BIN_OP_INT64)     return GAL_TYPE_INT64;
-          else if(GAL_CONFIG_BIN_OP_FLOAT32)   return GAL_TYPE_FLOAT32;
-          else if(GAL_CONFIG_BIN_OP_FLOAT64)   return GAL_TYPE_FLOAT64;
-        }
-      break;
 
-    case GAL_TYPE_UINT64:
-      if(GAL_CONFIG_BIN_OP_UINT64) return GAL_TYPE_UINT64;
-      else
-        {
-          if     (GAL_CONFIG_BIN_OP_FLOAT32)   return GAL_TYPE_FLOAT32;
-          else if(GAL_CONFIG_BIN_OP_FLOAT64)   return GAL_TYPE_FLOAT64;
-        }
-      break;
+  /* Clean up if necessary. Note that if the operation was requested to be
+     in place, then the output might be one of the inputs. */
+  if(flags & GAL_ARITHMETIC_FREE)
+    {
+      if     (o==l)       gal_data_free(r);
+      else if(o==r)       gal_data_free(l);
+      else              { gal_data_free(l); gal_data_free(r); }
+    }
 
-    case GAL_TYPE_INT64:
-      if(GAL_CONFIG_BIN_OP_INT64) return GAL_TYPE_INT64;
-      else
-        {
-          if     (GAL_CONFIG_BIN_OP_FLOAT32)   return GAL_TYPE_FLOAT32;
-          else if(GAL_CONFIG_BIN_OP_FLOAT64)   return GAL_TYPE_FLOAT64;
-        }
-      break;
 
-    case GAL_TYPE_FLOAT32:
-      if(GAL_CONFIG_BIN_OP_FLOAT32) return GAL_TYPE_FLOAT32;
-      else
-        {
-          if     (GAL_CONFIG_BIN_OP_FLOAT64)   return GAL_TYPE_FLOAT64;
-        }
-      break;
+  /* Return */
+  return o;
+}
 
-    case GAL_TYPE_FLOAT64:
-      if         (GAL_CONFIG_BIN_OP_FLOAT64)   return GAL_TYPE_FLOAT64;
-      break;
 
+
+
+
+#define BINFUNC_RUN_FUNCTION(OT, RT, LT, OP){                           \
+    LT *la=l->array;                                                    \
+    RT *ra=r->array;                                                    \
+    OT *oa=o->array, *of=oa + o->size;                                  \
+    if(l->size==r->size) do *oa = OP(*la++, *ra++); while(++oa<of);     \
+    else if(l->size==1)  do *oa = OP(*la,   *ra++); while(++oa<of);     \
+    else                 do *oa = OP(*la++, *ra  ); while(++oa<of);     \
+  }
+
+
+#define BINFUNC_F_OPERATOR_LEFT_RIGHT_SET(RT, LT, OP)                   \
+  switch(o->type)                                                       \
+    {                                                                   \
+    case GAL_TYPE_FLOAT32:                                              \
+      BINFUNC_RUN_FUNCTION(float, RT, LT, OP);                          \
+      break;                                                            \
+    case GAL_TYPE_FLOAT64:                                              \
+      BINFUNC_RUN_FUNCTION(double, RT, LT, OP);                         \
+      break;                                                            \
+    default:                                                            \
+      error(EXIT_FAILURE, 0, "%s: type %d not recognized for o->type ", \
+            "BINFUNC_F_OPERATOR_LEFT_RIGHT_SET", o->type);              \
+    }
+
+
+#define BINFUNC_F_OPERATOR_LEFT_SET(LT, OP)                             \
+  switch(r->type)                                                       \
+    {                                                                   \
+    case GAL_TYPE_FLOAT32:                                              \
+      BINFUNC_F_OPERATOR_LEFT_RIGHT_SET(float, LT, OP);                 \
+      break;                                                            \
+    case GAL_TYPE_FLOAT64:                                              \
+      BINFUNC_F_OPERATOR_LEFT_RIGHT_SET(double, LT, OP);                \
+      break;                                                            \
+    default:                                                            \
+      error(EXIT_FAILURE, 0, "%s: type %d not recognized for r->type",  \
+            "BINFUNC_F_OPERATOR_LEFT_SET", r->type);                    \
+    }
+
+
+#define BINFUNC_F_OPERATOR_SET(OP)                                      \
+  switch(l->type)                                                       \
+    {                                                                   \
+    case GAL_TYPE_FLOAT32:                                              \
+      BINFUNC_F_OPERATOR_LEFT_SET(float, OP);                           \
+      break;                                                            \
+    case GAL_TYPE_FLOAT64:                                              \
+      BINFUNC_F_OPERATOR_LEFT_SET(double, OP);                          \
+      break;                                                            \
+    default:                                                            \
+      error(EXIT_FAILURE, 0, "%s: type %d not recognized for l->type",  \
+            "BINFUNC_F_OPERATOR_SET", l->type);                         \
+    }
+
+
+static gal_data_t *
+arithmetic_binary_function_flt(int operator, int flags, gal_data_t *l,
+                               gal_data_t *r)
+{
+  int final_otype;
+  gal_data_t *o=NULL;
+  size_t out_size, minmapsize;
+
+
+  /* Simple sanity check on the input sizes */
+  if( !( (flags & GAL_ARITHMETIC_NUMOK) && (l->size==1 || r->size==1))
+      && gal_data_dsize_is_different(l, r) )
+    error(EXIT_FAILURE, 0, "%s: the input datasets don't have the same "
+          "dimension/size", __func__);
+
+  /* Check for the types of the left and right operands. */
+  arithmetic_check_float_input(l, operator, "first");
+  arithmetic_check_float_input(r, operator, "second");
+
+  /* Set the output type. */
+  final_otype = gal_type_out(l->type, r->type);
+
+  /* Set the output sizes. */
+  minmapsize = ( l->minmapsize < r->minmapsize
+                 ? l->minmapsize : r->minmapsize );
+  out_size = l->size > r->size ? l->size : r->size;
+
+
+  /* If we want inplace output, set the output pointer to one input. Note
+     that the output type can be different from both inputs.  */
+  if(flags & GAL_ARITHMETIC_INPLACE)
+    {
+      if     (l->type==final_otype && out_size==l->size)   o = l;
+      else if(r->type==final_otype && out_size==r->size)   o = r;
+    }
+
+
+  /* If the output pointer was not set for any reason, allocate it. For
+     `mmapsize', note that since its `size_t', it will always be
+     Positive. The `-1' that is recommended to give when you want the value
+     in RAM is actually the largest possible memory location. So we just
+     have to choose the smaller minmapsize of the two to decide if the
+     output array should be in RAM or not. */
+  if(o==NULL)
+    o = gal_data_alloc(NULL, final_otype,
+                       l->size>1 ? l->ndim  : r->ndim,
+                       l->size>1 ? l->dsize : r->dsize,
+                       l->size>1 ? l->wcs : r->wcs, 0, minmapsize,
+                       NULL, NULL, NULL);
+
+
+  /* Start setting the operator and operands. */
+  switch(operator)
+    {
+    case GAL_ARITHMETIC_OP_POW:  BINFUNC_F_OPERATOR_SET( pow  ); break;
     default:
-      error(EXIT_FAILURE, 0, "%s: type %d not recognized", __func__, intype);
+      error(EXIT_FAILURE, 0, "%s: operator code %d not recognized",
+            __func__, operator);
     }
 
-  return 0;
+
+  /* Clean up. Note that if the input arrays can be freed, and any of right
+     or left arrays needed conversion, `BINFUNC_CONVERT_TO_COMPILED_TYPE'
+     has already freed the input arrays, and we only have `r' and `l'
+     allocated in any case. Alternatively, when the inputs shouldn't be
+     freed, the only allocated spaces are the `r' and `l' arrays if their
+     types weren't compiled for binary operations, we can tell this from
+     the pointers: if they are different from the original pointers, they
+     were allocated. */
+  if(flags & GAL_ARITHMETIC_FREE)
+    {
+      if     (o==l)       gal_data_free(r);
+      else if(o==r)       gal_data_free(l);
+      else              { gal_data_free(l); gal_data_free(r); }
+    }
+
+  /* Return */
+  return o;
 }
 
 
@@ -1320,92 +1399,6 @@ gal_arithmetic_operator_string(int operator)
 
 
 
-/* Note that for signed types, we won't be considering the unsigned types
-   of the larger types. */
-gal_data_t *
-gal_arithmetic_convert_to_compiled_type(gal_data_t *in, int flags)
-{
-  int ntype;
-  char *typestring;
-  gal_data_t *out=NULL;
-
-  /* Set the best compiled type. */
-  ntype=arithmetic_nearest_compiled_type(in->type);
-
-  /* If type is not compiled, then convert the dataset to the first
-     compiled larger type. */
-  if(in->type==ntype)
-    out=in;
-  else
-    {
-      if(ntype)
-        {
-          out=gal_data_copy_to_new_type(in, ntype);
-          if(flags & GAL_ARITHMETIC_FREE)
-            { gal_data_free(in); in=NULL; }
-        }
-      else
-        {
-          typestring=gal_type_name(in->type, 1);
-          error(EXIT_FAILURE, 0, "The given %s type data given to "
-                "binary operators is not compiled for native operation "
-                "and no larger types are compiled either.\n\nThe "
-                "largest type (which can act as a fallback for any "
-                "input type is double, so configure Gnuastro again "
-                "with `--enable-bin-op-double' to not get this error "
-                "any more. However, if you commonly deal with %s type "
-                "data, also enable %s with a similar option at "
-                "configure time to greatly decrease running time and "
-                "avoid unnecessary RAM and CPU resources. Run"
-                "`./configure --help' in Gnuastro's top source "
-                "directory (after unpacking the tarball) for the full "
-                "list of options", typestring, typestring, typestring);
-        }
-    }
-
-  /* Return the output data structure */
-  if(out==NULL)
-    error(EXIT_FAILURE, 0, "%s: a bug! Please contact us at %s, so we can "
-          "fix the problem. For some reason, the `out' array has not been "
-          "set", __func__, PACKAGE_BUGREPORT);
-  return out;
-}
-
-
-
-
-
-/* Call functions in the `gnuastro/statistics' library. */
-static gal_data_t *
-arithmetic_from_statistics(int operator, int flags, gal_data_t *input)
-{
-  gal_data_t *out=NULL;
-  int ip=(flags & GAL_ARITHMETIC_INPLACE) || (flags & GAL_ARITHMETIC_FREE);
-
-  switch(operator)
-    {
-    case GAL_ARITHMETIC_OP_MINVAL:  out=gal_statistics_minimum(input); break;
-    case GAL_ARITHMETIC_OP_MAXVAL:  out=gal_statistics_maximum(input); break;
-    case GAL_ARITHMETIC_OP_NUMVAL:  out=gal_statistics_number(input);  break;
-    case GAL_ARITHMETIC_OP_SUMVAL:  out=gal_statistics_sum(input);     break;
-    case GAL_ARITHMETIC_OP_MEANVAL: out=gal_statistics_mean(input);    break;
-    case GAL_ARITHMETIC_OP_STDVAL:  out=gal_statistics_std(input);     break;
-    case GAL_ARITHMETIC_OP_MEDIANVAL:
-      out=gal_statistics_median(input, ip); break;
-    default:
-      error(EXIT_FAILURE, 0, "%s: operator code %d not recognized",
-            __func__, operator);
-    }
-
-  /* If the input is to be freed, then do so and return the output. */
-  if( flags & GAL_ARITHMETIC_FREE ) gal_data_free(input);
-  return out;
-}
-
-
-
-
-
 gal_data_t *
 gal_arithmetic(int operator, int flags, ...)
 {
@@ -1514,12 +1507,12 @@ gal_arithmetic(int operator, int flags, ...)
     case GAL_ARITHMETIC_OP_MODULO:
       d1 = va_arg(va, gal_data_t *);
       d2 = va_arg(va, gal_data_t *);
-      out=arithmetic_onlyint_binary(operator, flags, d1, d2);
+      out=arithmetic_binary(operator, flags, d1, d2);
       break;
 
     case GAL_ARITHMETIC_OP_BITNOT:
       d1 = va_arg(va, gal_data_t *);
-      out=arithmetic_onlyint_bitwise_not(flags, d1);
+      out=arithmetic_bitwise_not(flags, d1);
       break;
 
 
diff --git a/lib/gnuastro-internal/arithmetic-onlyint.h 
b/lib/gnuastro-internal/arithmetic-and.h
similarity index 78%
rename from lib/gnuastro-internal/arithmetic-onlyint.h
rename to lib/gnuastro-internal/arithmetic-and.h
index 2a8b3bd..4632a16 100644
--- a/lib/gnuastro-internal/arithmetic-onlyint.h
+++ b/lib/gnuastro-internal/arithmetic-and.h
@@ -20,16 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __GAL_ARITHMETIC_ONLYINT_H__
-#define __GAL_ARITHMETIC_ONLYINT_H__
-
-
-gal_data_t *
-arithmetic_onlyint_binary(int operator, int flags, gal_data_t *lo,
-                          gal_data_t *ro);
-
-gal_data_t *
-arithmetic_onlyint_bitwise_not(int flags, gal_data_t *in);
+#ifndef __ARITHMETIC_AND_H__
+#define __ARITHMETIC_AND_H__
 
+void
+arithmetic_and(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-binary.h
index c441a4a..b1b79c3 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-binary.h
@@ -24,8 +24,254 @@ along with Gnuastro. If not, see 
<http://www.gnu.org/licenses/>.
 #define __ARITHMETIC_BINARY_H__
 
 
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+
+
+
+
+
+
+
+
+/************************************************************************/
+/*************             Low-level operators          *****************/
+/************************************************************************/
+/* Final step to be used by all operators and all types. */
+#define BINARY_OP_OT_RT_LT_SET(OP, OT, LT, RT) {                        \
+    LT lb, *la=l->array;                                                \
+    RT rb, *ra=r->array;                                                \
+    OT ob, *oa=o->array, *of=oa + o->size;                              \
+    if(checkblank)                                                      \
+      {                                                                 \
+        gal_blank_write(&lb, l->type);                                  \
+        gal_blank_write(&rb, r->type);                                  \
+        gal_blank_write(&ob, o->type);                                  \
+        do                                                              \
+          {                                                             \
+            if(lb==lb && rb==rb)/* Both are integers.                */ \
+              *oa = (*la!=lb  && *ra!=rb)  ? *la OP *ra : ob ;          \
+            else if(lb==lb)     /* Only left operand is an integer.  */ \
+              *oa = (*la!=lb  && *ra==*ra) ? *la OP *ra : ob;           \
+            else                /* Only right operand is an integer. */ \
+              *oa = (*la==*la && *ra!=rb)  ? *la OP *ra : ob;           \
+            if(l->size>1) ++la;                                         \
+            if(r->size>1) ++ra;                                         \
+          }                                                             \
+        while(++oa<of);                                                 \
+      }                                                                 \
+    else                                                                \
+      {                                                                 \
+        if(l->size==r->size) do *oa = *la++ OP *ra++; while(++oa<of);   \
+        else if(l->size==1)  do *oa = *la   OP *ra++; while(++oa<of);   \
+        else                 do *oa = *la++ OP *ra;   while(++oa<of);   \
+      }                                                                 \
+  }
+
+
+
+
+
+/* Blank values aren't defined for integer operators. */
+#define BINARY_INT_OP_OT_RT_LT_SET(OP, OT, LT, RT) {               \
+    LT *la=l->array;                                               \
+    RT *ra=r->array;                                               \
+    OT *oa=o->array, *of=oa + o->size;                             \
+    if(l->size==r->size) do *oa = *la++ OP *ra++; while(++oa<of);  \
+    else if(l->size==1)  do *oa = *la   OP *ra++; while(++oa<of);  \
+    else                 do *oa = *la++ OP *ra;   while(++oa<of);  \
+  }
+
+
+
+
+
+/* This is for operators like `&&' and `||', where the right operator is
+   not necessarily read (and thus incremented). */
+#define BINARY_OP_INCR_OT_RT_LT_SET(OP, OT, LT, RT) {                   \
+    LT *la=l->array;                                                    \
+    RT *ra=r->array;                                                    \
+    OT *oa=o->array, *of=oa + o->size;                                  \
+    if(l->size==r->size) do {*oa = *la++ OP *ra; ++ra;} while(++oa<of); \
+    else if(l->size==1)  do {*oa = *la   OP *ra; ++ra;} while(++oa<of); \
+    else                 do  *oa = *la++ OP *ra;        while(++oa<of); \
+  }
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+/************************************************************************/
+/*************              Type specifiers             *****************/
+/************************************************************************/
+
+/* Flags for BINARY_OP_RT_LT_SET to identify the output type. */
+enum arithmetic_binary_outtype_flags
+{
+  ARITHMETIC_BINARY_INVALID,                /* ==0 by C standard. */
+
+  ARITHMETIC_BINARY_OUT_TYPE_LEFT,
+  ARITHMETIC_BINARY_OUT_TYPE_RIGHT,
+  ARITHMETIC_BINARY_OUT_TYPE_UINT8,
+  ARITHMETIC_BINARY_OUT_TYPE_INCR_SEP,
+};
+
+
+
+
+
+/* For operators whose type may be any of the given inputs only for
+   integers. For integer operators we have less options. */
+#define BINARY_SET_OUT_INT(F, OP, LT, RT)                               \
+  switch(F)                                                             \
+    {                                                                   \
+    case ARITHMETIC_BINARY_OUT_TYPE_LEFT:                               \
+      BINARY_INT_OP_OT_RT_LT_SET(OP, LT, LT, RT);                       \
+      break;                                                            \
+    case ARITHMETIC_BINARY_OUT_TYPE_RIGHT:                              \
+      BINARY_INT_OP_OT_RT_LT_SET(OP, RT, LT, RT);                       \
+      break;                                                            \
+    default:                                                            \
+      error(EXIT_FAILURE, 0, "%s: a bug! please contact us at %s to "   \
+            "address the problem. %d not recognized for `F'",           \
+            "BINARY_SET_OUT_INT", PACKAGE_BUGREPORT, F);                \
+    }
+
+
+
+
+
+/* For operators whose type may be any of the given inputs. */
+#define BINARY_SET_OUT(F, OP, LT, RT)                                   \
+  switch(F)                                                             \
+    {                                                                   \
+    case ARITHMETIC_BINARY_OUT_TYPE_LEFT:                               \
+      BINARY_OP_OT_RT_LT_SET(OP,      LT,      LT, RT);                 \
+      break;                                                            \
+    case ARITHMETIC_BINARY_OUT_TYPE_RIGHT:                              \
+      BINARY_OP_OT_RT_LT_SET(OP,      RT,      LT, RT);                 \
+      break;                                                            \
+    case ARITHMETIC_BINARY_OUT_TYPE_UINT8:                              \
+      BINARY_OP_OT_RT_LT_SET(OP,      uint8_t, LT, RT);                 \
+      break;                                                            \
+    case ARITHMETIC_BINARY_OUT_TYPE_INCR_SEP:                           \
+      BINARY_OP_INCR_OT_RT_LT_SET(OP, uint8_t, LT, RT);                 \
+      break;                                                            \
+    default:                                                            \
+      error(EXIT_FAILURE, 0, "%s: a bug! please contact us at %s to "   \
+            "address the problem. %d not recognized for `F'",           \
+            "BINARY_SET_OUT", PACKAGE_BUGREPORT, F);                    \
+    }
+
+
+
+
+
+/* Set the right operator type only for integers (integer operators can't
+   take floating point types). So floating point types must not be in the
+   list of possibilities. */
+#define BINARY_SET_RT_INT(F, OP, LT)                                        \
+  switch(r->type)                                                           \
+    {                                                                       \
+    case GAL_TYPE_UINT8:   BINARY_SET_OUT_INT( F, OP, LT, uint8_t  ) break; \
+    case GAL_TYPE_INT8:    BINARY_SET_OUT_INT( F, OP, LT, int8_t   ) break; \
+    case GAL_TYPE_UINT16:  BINARY_SET_OUT_INT( F, OP, LT, uint16_t ) break; \
+    case GAL_TYPE_INT16:   BINARY_SET_OUT_INT( F, OP, LT, int16_t  ) break; \
+    case GAL_TYPE_UINT32:  BINARY_SET_OUT_INT( F, OP, LT, uint32_t ) break; \
+    case GAL_TYPE_INT32:   BINARY_SET_OUT_INT( F, OP, LT, int32_t  ) break; \
+    case GAL_TYPE_UINT64:  BINARY_SET_OUT_INT( F, OP, LT, uint64_t ) break; \
+    case GAL_TYPE_INT64:   BINARY_SET_OUT_INT( F, OP, LT, int64_t  ) break; \
+    default:                                                                \
+      error(EXIT_FAILURE, 0, "%s: a bug! Please contact us at %s to "       \
+            "address the problem. %d is not a usable type code",            \
+            "BINARY_SET_RT", PACKAGE_BUGREPORT, r->type);                   \
+    }
+
+
+
+
+
+/* Set the right operator type. */
+#define BINARY_SET_RT(F, OP, LT)                                        \
+  switch(r->type)                                                       \
+    {                                                                   \
+    case GAL_TYPE_UINT8:   BINARY_SET_OUT( F, OP, LT, uint8_t  ) break; \
+    case GAL_TYPE_INT8:    BINARY_SET_OUT( F, OP, LT, int8_t   ) break; \
+    case GAL_TYPE_UINT16:  BINARY_SET_OUT( F, OP, LT, uint16_t ) break; \
+    case GAL_TYPE_INT16:   BINARY_SET_OUT( F, OP, LT, int16_t  ) break; \
+    case GAL_TYPE_UINT32:  BINARY_SET_OUT( F, OP, LT, uint32_t ) break; \
+    case GAL_TYPE_INT32:   BINARY_SET_OUT( F, OP, LT, int32_t  ) break; \
+    case GAL_TYPE_UINT64:  BINARY_SET_OUT( F, OP, LT, uint64_t ) break; \
+    case GAL_TYPE_INT64:   BINARY_SET_OUT( F, OP, LT, int64_t  ) break; \
+    case GAL_TYPE_FLOAT32: BINARY_SET_OUT( F, OP, LT, float    ) break; \
+    case GAL_TYPE_FLOAT64: BINARY_SET_OUT( F, OP, LT, double   ) break; \
+    default:                                                            \
+      error(EXIT_FAILURE, 0, "%s: a bug! Please contact us at %s to "   \
+            "address the problem. %d is not a usable type code",        \
+            "BINARY_SET_RT", PACKAGE_BUGREPORT, r->type);               \
+    }
+
+
+
+
+
+/* Set the left operator type only for integers (integer operators can't
+   take floating point types). So floating point types must not be in the
+   list of possibilities. */
+#define BINARY_SET_LT_INT(F, OP)                                        \
+  switch(l->type)                                                       \
+    {                                                                   \
+    case GAL_TYPE_UINT8:   BINARY_SET_RT_INT( F, OP, uint8_t  ) break;  \
+    case GAL_TYPE_INT8:    BINARY_SET_RT_INT( F, OP, int8_t   ) break;  \
+    case GAL_TYPE_UINT16:  BINARY_SET_RT_INT( F, OP, uint16_t ) break;  \
+    case GAL_TYPE_INT16:   BINARY_SET_RT_INT( F, OP, int16_t  ) break;  \
+    case GAL_TYPE_UINT32:  BINARY_SET_RT_INT( F, OP, uint32_t ) break;  \
+    case GAL_TYPE_INT32:   BINARY_SET_RT_INT( F, OP, int32_t  ) break;  \
+    case GAL_TYPE_UINT64:  BINARY_SET_RT_INT( F, OP, uint64_t ) break;  \
+    case GAL_TYPE_INT64:   BINARY_SET_RT_INT( F, OP, int64_t  ) break;  \
+    default:                                                            \
+      error(EXIT_FAILURE, 0, "%s: a bug! Please contact us at %s to "   \
+            "address the problem. %d is not a usable type code",        \
+            "BINARY_SET_LT_INT", PACKAGE_BUGREPORT, l->type);           \
+    }
+
+
+
+
+
+/* Set the left operator type. */
+#define BINARY_SET_LT(F, OP)                                            \
+  switch(l->type)                                                       \
+    {                                                                   \
+    case GAL_TYPE_UINT8:   BINARY_SET_RT( F, OP, uint8_t  ) break;      \
+    case GAL_TYPE_INT8:    BINARY_SET_RT( F, OP, int8_t   ) break;      \
+    case GAL_TYPE_UINT16:  BINARY_SET_RT( F, OP, uint16_t ) break;      \
+    case GAL_TYPE_INT16:   BINARY_SET_RT( F, OP, int16_t  ) break;      \
+    case GAL_TYPE_UINT32:  BINARY_SET_RT( F, OP, uint32_t ) break;      \
+    case GAL_TYPE_INT32:   BINARY_SET_RT( F, OP, int32_t  ) break;      \
+    case GAL_TYPE_UINT64:  BINARY_SET_RT( F, OP, uint64_t ) break;      \
+    case GAL_TYPE_INT64:   BINARY_SET_RT( F, OP, int64_t  ) break;      \
+    case GAL_TYPE_FLOAT32: BINARY_SET_RT( F, OP, float    ) break;      \
+    case GAL_TYPE_FLOAT64: BINARY_SET_RT( F, OP, double   ) break;      \
+    default:                                                            \
+      error(EXIT_FAILURE, 0, "%s: a bug! Please contact us at %s to "   \
+            "address the problem. %d is not a usable type code",        \
+            "BINARY_SET_LT", PACKAGE_BUGREPORT, l->type);               \
+    }
+
 
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-bitand.h
similarity index 86%
copy from lib/gnuastro-internal/arithmetic-binary.h
copy to lib/gnuastro-internal/arithmetic-bitand.h
index c441a4a..bb221a0 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-bitand.h
@@ -20,12 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __ARITHMETIC_BINARY_H__
-#define __ARITHMETIC_BINARY_H__
-
-
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+#ifndef __ARITHMETIC_BITAND_H__
+#define __ARITHMETIC_BITAND_H__
 
+void
+arithmetic_bitand(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-bitlsh.h
similarity index 86%
copy from lib/gnuastro-internal/arithmetic-binary.h
copy to lib/gnuastro-internal/arithmetic-bitlsh.h
index c441a4a..8311124 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-bitlsh.h
@@ -20,12 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __ARITHMETIC_BINARY_H__
-#define __ARITHMETIC_BINARY_H__
-
-
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+#ifndef __ARITHMETIC_BITLSH_H__
+#define __ARITHMETIC_BITLSH_H__
 
+void
+arithmetic_bitlsh(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-bitor.h
similarity index 86%
copy from lib/gnuastro-internal/arithmetic-binary.h
copy to lib/gnuastro-internal/arithmetic-bitor.h
index c441a4a..5cd2cb5 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-bitor.h
@@ -20,12 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __ARITHMETIC_BINARY_H__
-#define __ARITHMETIC_BINARY_H__
-
-
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+#ifndef __ARITHMETIC_BITOR_H__
+#define __ARITHMETIC_BITOR_H__
 
+void
+arithmetic_bitor(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-bitrsh.h
similarity index 86%
copy from lib/gnuastro-internal/arithmetic-binary.h
copy to lib/gnuastro-internal/arithmetic-bitrsh.h
index c441a4a..bbe8695 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-bitrsh.h
@@ -20,12 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __ARITHMETIC_BINARY_H__
-#define __ARITHMETIC_BINARY_H__
-
-
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+#ifndef __ARITHMETIC_BITRSH_H__
+#define __ARITHMETIC_BITRSH_H__
 
+void
+arithmetic_bitrsh(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-bitxor.h
similarity index 86%
copy from lib/gnuastro-internal/arithmetic-binary.h
copy to lib/gnuastro-internal/arithmetic-bitxor.h
index c441a4a..a94eede 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-bitxor.h
@@ -20,12 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __ARITHMETIC_BINARY_H__
-#define __ARITHMETIC_BINARY_H__
-
-
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+#ifndef __ARITHMETIC_BITXOR_H__
+#define __ARITHMETIC_BITXOR_H__
 
+void
+arithmetic_bitxor(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-divide.h
similarity index 86%
copy from lib/gnuastro-internal/arithmetic-binary.h
copy to lib/gnuastro-internal/arithmetic-divide.h
index c441a4a..fb93428 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-divide.h
@@ -20,12 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __ARITHMETIC_BINARY_H__
-#define __ARITHMETIC_BINARY_H__
-
-
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+#ifndef __ARITHMETIC_DIVIDE_H__
+#define __ARITHMETIC_DIVIDE_H__
 
+void
+arithmetic_divide(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-eq.h
similarity index 86%
copy from lib/gnuastro-internal/arithmetic-binary.h
copy to lib/gnuastro-internal/arithmetic-eq.h
index c441a4a..a5947bf 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-eq.h
@@ -20,12 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __ARITHMETIC_BINARY_H__
-#define __ARITHMETIC_BINARY_H__
-
-
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+#ifndef __ARITHMETIC_EQ_H__
+#define __ARITHMETIC_EQ_H__
 
+void
+arithmetic_eq(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-ge.h
similarity index 86%
copy from lib/gnuastro-internal/arithmetic-binary.h
copy to lib/gnuastro-internal/arithmetic-ge.h
index c441a4a..f1346f7 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-ge.h
@@ -20,12 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __ARITHMETIC_BINARY_H__
-#define __ARITHMETIC_BINARY_H__
-
-
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+#ifndef __ARITHMETIC_GE_H__
+#define __ARITHMETIC_GE_H__
 
+void
+arithmetic_ge(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-gt.h
similarity index 86%
copy from lib/gnuastro-internal/arithmetic-binary.h
copy to lib/gnuastro-internal/arithmetic-gt.h
index c441a4a..5b2c9e8 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-gt.h
@@ -20,12 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __ARITHMETIC_BINARY_H__
-#define __ARITHMETIC_BINARY_H__
-
-
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+#ifndef __ARITHMETIC_GT_H__
+#define __ARITHMETIC_GT_H__
 
+void
+arithmetic_gt(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-internal.h 
b/lib/gnuastro-internal/arithmetic-internal.h
index 58c9041..d9cd846 100644
--- a/lib/gnuastro-internal/arithmetic-internal.h
+++ b/lib/gnuastro-internal/arithmetic-internal.h
@@ -59,7 +59,7 @@ along with Gnuastro. If not, see 
<http://www.gnu.org/licenses/>.
 __BEGIN_C_DECLS  /* From C++ preparations */
 
 int
-gal_arithmetic_binary_out_type(int operator, gal_data_t *l, gal_data_t *r);
+gal_arithmetic_binary_checkblank(gal_data_t *l, gal_data_t *r);
 
 char *
 gal_arithmetic_operator_string(int operator);
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-le.h
similarity index 86%
copy from lib/gnuastro-internal/arithmetic-binary.h
copy to lib/gnuastro-internal/arithmetic-le.h
index c441a4a..b45866e 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-le.h
@@ -20,12 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __ARITHMETIC_BINARY_H__
-#define __ARITHMETIC_BINARY_H__
-
-
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+#ifndef __ARITHMETIC_LE_H__
+#define __ARITHMETIC_LE_H__
 
+void
+arithmetic_le(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-lt.h
similarity index 86%
copy from lib/gnuastro-internal/arithmetic-binary.h
copy to lib/gnuastro-internal/arithmetic-lt.h
index c441a4a..40fe879 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-lt.h
@@ -20,12 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __ARITHMETIC_BINARY_H__
-#define __ARITHMETIC_BINARY_H__
-
-
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+#ifndef __ARITHMETIC_LT_H__
+#define __ARITHMETIC_LT_H__
 
+void
+arithmetic_lt(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-minus.h
similarity index 86%
copy from lib/gnuastro-internal/arithmetic-binary.h
copy to lib/gnuastro-internal/arithmetic-minus.h
index c441a4a..5a977b0 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-minus.h
@@ -20,12 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __ARITHMETIC_BINARY_H__
-#define __ARITHMETIC_BINARY_H__
-
-
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+#ifndef __ARITHMETIC_MINUS_H__
+#define __ARITHMETIC_MINUS_H__
 
+void
+arithmetic_minus(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-modulo.h
similarity index 86%
copy from lib/gnuastro-internal/arithmetic-binary.h
copy to lib/gnuastro-internal/arithmetic-modulo.h
index c441a4a..2671d5f 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-modulo.h
@@ -20,12 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __ARITHMETIC_BINARY_H__
-#define __ARITHMETIC_BINARY_H__
-
-
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+#ifndef __ARITHMETIC_MODULO_H__
+#define __ARITHMETIC_MODULO_H__
 
+void
+arithmetic_modulo(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-multiply.h
similarity index 86%
copy from lib/gnuastro-internal/arithmetic-binary.h
copy to lib/gnuastro-internal/arithmetic-multiply.h
index c441a4a..b951006 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-multiply.h
@@ -20,12 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __ARITHMETIC_BINARY_H__
-#define __ARITHMETIC_BINARY_H__
-
-
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+#ifndef __ARITHMETIC_MULTIPLY_H__
+#define __ARITHMETIC_MULTIPLY_H__
 
+void
+arithmetic_multiply(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-ne.h
similarity index 86%
copy from lib/gnuastro-internal/arithmetic-binary.h
copy to lib/gnuastro-internal/arithmetic-ne.h
index c441a4a..e54b9c0 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-ne.h
@@ -20,12 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __ARITHMETIC_BINARY_H__
-#define __ARITHMETIC_BINARY_H__
-
-
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+#ifndef __ARITHMETIC_NE_H__
+#define __ARITHMETIC_NE_H__
 
+void
+arithmetic_ne(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-or.h
similarity index 86%
copy from lib/gnuastro-internal/arithmetic-binary.h
copy to lib/gnuastro-internal/arithmetic-or.h
index c441a4a..e2763f6 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-or.h
@@ -20,12 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __ARITHMETIC_BINARY_H__
-#define __ARITHMETIC_BINARY_H__
-
-
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+#ifndef __ARITHMETIC_OR_H__
+#define __ARITHMETIC_OR_H__
 
+void
+arithmetic_or(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif
diff --git a/lib/gnuastro-internal/arithmetic-binary.h 
b/lib/gnuastro-internal/arithmetic-plus.h
similarity index 86%
copy from lib/gnuastro-internal/arithmetic-binary.h
copy to lib/gnuastro-internal/arithmetic-plus.h
index c441a4a..cf6ddcc 100644
--- a/lib/gnuastro-internal/arithmetic-binary.h
+++ b/lib/gnuastro-internal/arithmetic-plus.h
@@ -20,12 +20,10 @@ General Public License for more details.
 You should have received a copy of the GNU General Public License
 along with Gnuastro. If not, see <http://www.gnu.org/licenses/>.
 **********************************************************************/
-#ifndef __ARITHMETIC_BINARY_H__
-#define __ARITHMETIC_BINARY_H__
-
-
-gal_data_t *
-arithmetic_binary(int operator, int flags, gal_data_t *lo, gal_data_t *ro);
+#ifndef __ARITHMETIC_PLUS_H__
+#define __ARITHMETIC_PLUS_H__
 
+void
+arithmetic_plus(gal_data_t *l, gal_data_t *r, gal_data_t *o);
 
 #endif



reply via email to

[Prev in Thread] Current Thread [Next in Thread]