Support for OpenMP added

The configure script is now able to detect OpenMP and use it for
compilation and linking. The function hpix_bmp_projection_trace has
been modified to use OpenMP: the speed up is not large, but it is a
nice way to check that OpenMP works within HPixLib.

.gitignore

@@ -11,6 +11,6 @@ configure
 depcomp
 install-sh
 ltmain.sh
-m4/
+m4/l*
 missing
 src/config.h.in

configure.ac

@@ -20,8 +20,45 @@ PKG_CHECK_MODULES([cairo], [cairo],
 
 AC_CHECK_LIB(cfitsio, ffopen,, AC_MSG_ERROR(Cannot find the CFITSIO library.))
 
+########################################################################
+# Check for OpenMP
+#
+# The logic of the code is as follows:
+#
+# 1. If the user specified --with-openmp, force a check and produce
+#    an error message if OpenMP is not available
+#
+# 2. If the user did not specify --with-openmp, try to use OpenMP
+#    but only print a warning if it is not available
+
+AC_ARG_WITH([openmp],
+	[AS_HELP_STRING([--with-openmp], 
+			[Multi-threaded parallel version using OpenMP])],
+	[],
+	[with_openmp=check])
+
+AS_IF([test "x$with_openmp" != xno],
+      [AX_OPENMP([openmp=yes
+                  AC_DEFINE(HAVE_OPENMP,
+                            1,
+                            [Define to 1 if you have OpenMP available])
+                  CFLAGS="$CFLAGS $OPENMP_CFLAGS"
+                 ],
+                 [openmp=no
+                  AC_MSG_WARN([The C compiler does not support OpenMP])
+                  if test "x$with_openmp" != xcheck; then
+                      AC_MSG_FAILURE([--with-openmp was specified, but no OpenMP support is available])
+                  fi
+                 ])
+      ],
+      [openmp=no])
+
+########################################################################
+
 AM_CONDITIONAL([CAIRO_PRESENT], [test x$cairo = xyes])
+AM_CONDITIONAL([OPENMP_PRESENT], [test x$openmp = xyes])
 AC_SUBST(havecairo, $cairo)
+AC_SUBST(haveopenmp, $openmp)
 
 AC_OUTPUT(Makefile
 	src/Makefile
@@ -35,4 +72,5 @@ echo ""
 echo "  HPixLib will be compiled with the following options:"
 echo ""
 echo "  cairo: $cairo"
+echo "  OpenMP: $openmp"
 echo ""

examples/Makefile.am

@@ -30,14 +30,12 @@
 # as representing official policies, either expressed or implied, of
 # Maurizio Tomasi.
 
-AM_CFLAGS = -std=c99
-
 bin_PROGRAMS = map2ppm mapinfo
 
 AM_CPPFLAGS = -I$(top_srcdir)/src
 
 AM_CFLAGS += $(cairo_CFLAGS)
-LIBS += $(cairo_LIBS) -lcfitsio
+LIBS += $(cairo_LIBS) -lcfitsio -lm
 
 map2ppm_SOURCES = map2ppm.c
 map2ppm_LDADD = ../src/libhpix.la

m4/ax_openmp.m4

@@ -0,0 +1,109 @@
+# ===========================================================================
+#         http://www.gnu.org/software/autoconf-archive/ax_openmp.html
+# ===========================================================================
+#
+# SYNOPSIS
+#
+#   AX_OPENMP([ACTION-IF-FOUND[, ACTION-IF-NOT-FOUND]])
+#
+# DESCRIPTION
+#
+#   This macro tries to find out how to compile programs that use OpenMP a
+#   standard API and set of compiler directives for parallel programming
+#   (see http://www-unix.mcs/)
+#
+#   On success, it sets the OPENMP_CFLAGS/OPENMP_CXXFLAGS/OPENMP_F77FLAGS
+#   output variable to the flag (e.g. -omp) used both to compile *and* link
+#   OpenMP programs in the current language.
+#
+#   NOTE: You are assumed to not only compile your program with these flags,
+#   but also link it with them as well.
+#
+#   If you want to compile everything with OpenMP, you should set:
+#
+#     CFLAGS="$CFLAGS $OPENMP_CFLAGS"
+#     #OR#  CXXFLAGS="$CXXFLAGS $OPENMP_CXXFLAGS"
+#     #OR#  FFLAGS="$FFLAGS $OPENMP_FFLAGS"
+#
+#   (depending on the selected language).
+#
+#   The user can override the default choice by setting the corresponding
+#   environment variable (e.g. OPENMP_CFLAGS).
+#
+#   ACTION-IF-FOUND is a list of shell commands to run if an OpenMP flag is
+#   found, and ACTION-IF-NOT-FOUND is a list of commands to run it if it is
+#   not found. If ACTION-IF-FOUND is not specified, the default action will
+#   define HAVE_OPENMP.
+#
+# LICENSE
+#
+#   Copyright (c) 2008 Steven G. Johnson <stevenj@alum.mit.edu>
+#
+#   This program 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.
+#
+#   This program 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 this program. If not, see <http://www.gnu.org/licenses/>.
+#
+#   As a special exception, the respective Autoconf Macro's copyright owner
+#   gives unlimited permission to copy, distribute and modify the configure
+#   scripts that are the output of Autoconf when processing the Macro. You
+#   need not follow the terms of the GNU General Public License when using
+#   or distributing such scripts, even though portions of the text of the
+#   Macro appear in them. The GNU General Public License (GPL) does govern
+#   all other use of the material that constitutes the Autoconf Macro.
+#
+#   This special exception to the GPL applies to versions of the Autoconf
+#   Macro released by the Autoconf Archive. When you make and distribute a
+#   modified version of the Autoconf Macro, you may extend this special
+#   exception to the GPL to apply to your modified version as well.
+
+#serial 9
+
+AC_DEFUN([AX_OPENMP], [
+AC_PREREQ(2.59) dnl for _AC_LANG_PREFIX
+
+AC_CACHE_CHECK([for OpenMP flag of _AC_LANG compiler], ax_cv_[]_AC_LANG_ABBREV[]_openmp, [save[]_AC_LANG_PREFIX[]FLAGS=$[]_AC_LANG_PREFIX[]FLAGS
+ax_cv_[]_AC_LANG_ABBREV[]_openmp=unknown
+# Flags to try:  -fopenmp (gcc), -openmp (icc), -mp (SGI & PGI),
+#                -xopenmp (Sun), -omp (Tru64), -qsmp=omp (AIX), none
+ax_openmp_flags="-fopenmp -openmp -mp -xopenmp -omp -qsmp=omp none"
+if test "x$OPENMP_[]_AC_LANG_PREFIX[]FLAGS" != x; then
+  ax_openmp_flags="$OPENMP_[]_AC_LANG_PREFIX[]FLAGS $ax_openmp_flags"
+fi
+for ax_openmp_flag in $ax_openmp_flags; do
+  case $ax_openmp_flag in
+    none) []_AC_LANG_PREFIX[]FLAGS=$save[]_AC_LANG_PREFIX[] ;;
+    *) []_AC_LANG_PREFIX[]FLAGS="$save[]_AC_LANG_PREFIX[]FLAGS $ax_openmp_flag" ;;
+  esac
+  AC_TRY_LINK([#ifdef __cplusplus
+extern "C"
+#endif
+void omp_set_num_threads(int);], [const int N = 100000;
+  int i, arr[N];
+
+  omp_set_num_threads(2);
+
+  #pragma omp parallel for
+  for (i = 0; i < N; i++) {
+    arr[i] = i;
+  }], [ax_cv_[]_AC_LANG_ABBREV[]_openmp=$ax_openmp_flag; break])
+done
+[]_AC_LANG_PREFIX[]FLAGS=$save[]_AC_LANG_PREFIX[]FLAGS
+])
+if test "x$ax_cv_[]_AC_LANG_ABBREV[]_openmp" = "xunknown"; then
+  m4_default([$2],:)
+else
+  if test "x$ax_cv_[]_AC_LANG_ABBREV[]_openmp" != "xnone"; then
+    OPENMP_[]_AC_LANG_PREFIX[]FLAGS=$ax_cv_[]_AC_LANG_ABBREV[]_openmp
+  fi
+  m4_default([$1], [AC_DEFINE(HAVE_OPENMP,1,[Define if OpenMP is enabled])])
+fi
+])dnl AX_OPENMP

src/Makefile.am

@@ -30,8 +30,6 @@
 # as representing official policies, either expressed or implied, of
 # Maurizio Tomasi.
 
-AM_CFLAGS = -std=c99
-
 LIBRARIES_TO_BUILD = libhpix.la
 if CAIRO_PRESENT
 LIBRARIES_TO_BUILD += libhpix_cairo.la

src/bitmap.c

@@ -18,6 +18,7 @@
 #include "config.h"
 
 #include <hpixlib/hpix.h>
+#include <float.h>
 #include <math.h>
 #include <assert.h>
 
@@ -175,49 +176,56 @@ hpix_bmp_projection_trace(const hpix_bmp_projection_t * proj,
 	? hpix_angles_to_nest_pixel
 	: hpix_angles_to_ring_pixel;
 
+    size_t num_of_pixels = proj->width * proj->height;
     double *restrict bitmap =
-	hpix_malloc(sizeof(bitmap[0]), proj->width * proj->height);
-    /* Although the two pointers point to the same area, only one of
-     * them (bitmap_ptr) will be used for writing. So we're justified
-     * in using "restricted" for them. */
-    double *restrict bitmap_ptr = bitmap;
+	hpix_malloc(sizeof(bitmap[0]), num_of_pixels);
 
     double * pixels = hpix_map_pixels(map);
-    int minmax_flag = 0;
+
+    /* First step: render the bitmap */
+#pragma omp parallel for default(shared)
     for (unsigned int y = 0; y < hpix_bmp_projection_height(proj); ++y)
     {
-	for (unsigned int x = 0; x < hpix_bmp_projection_width(proj); ++x)
+	double * line_ptr = bitmap + y * hpix_bmp_projection_width(proj);
+
+	for (unsigned int x = 0; 
+	     x < hpix_bmp_projection_width(proj); 
+	     ++x, ++line_ptr)
 	{
 	    double theta, phi;
 
 	    if(! hpix_bmp_projection_xy_to_angles(proj, x, y, &theta, &phi))
 	    {
-		*bitmap_ptr++ = INFINITY; /* Skip the pixel */
+		*line_ptr = INFINITY; /* Skip the pixel */
 		continue;
 	    }
 
 	    hpix_pixel_num_t pixel_idx =
 		angles_to_pixel_fn(nside, theta, phi);
 	    if(pixels[pixel_idx] > -1.6e+30)
-		*bitmap_ptr = pixels[pixel_idx];
+		*line_ptr = pixels[pixel_idx];
 	    else
-		*bitmap_ptr = NAN;
+		*line_ptr = NAN;
+	}
+    }
 
-	    if(! minmax_flag)
-	    {
-		minmax_flag = 1;
-		if(min_value)
-		    *min_value = *bitmap_ptr;
-		if(max_value)
-		    *max_value = *bitmap_ptr;
-	    } else {
-		if(min_value && *min_value > *bitmap_ptr)
-		    *min_value = *bitmap_ptr;
-		if(max_value && *max_value < *bitmap_ptr)
-		    *max_value = *bitmap_ptr;
-	    }
+    /* Second step: if the user asked for them, compute the maximum
+     * and/or minimum values in the bitmap */
+    if(min_value || max_value) 
+    {
+	if(min_value)
+	    *min_value = FLT_MAX;
 
-	    bitmap_ptr++;
+	if(max_value)
+	    *max_value = FLT_MIN;
+
+	double * bitmap_ptr = bitmap;
+	for(size_t idx = 0; idx < num_of_pixels; ++idx, ++bitmap_ptr) 
+	{
+	    if(min_value && *min_value > *bitmap_ptr)
+		*min_value = *bitmap_ptr;
+	    if(max_value && *max_value < *bitmap_ptr)
+		*max_value = *bitmap_ptr;
 	}
     }
 

utilities/Makefile.am

@@ -30,8 +30,6 @@
 # as representing official policies, either expressed or implied, of
 # Maurizio Tomasi.
 
-AM_CFLAGS = -std=c99
-
 if CAIRO_PRESENT
 PROGRAMS_TO_BUILD = map2fig
 else
@@ -42,7 +40,7 @@ bin_PROGRAMS = $(PROGRAMS_TO_BUILD)
 
 AM_CPPFLAGS = -I$(top_srcdir)/src
 
-AM_CFLAGS += $(cairo_CFLAGS)
+AM_CFLAGS = $(cairo_CFLAGS)
 LIBS += $(cairo_LIBS) -lcfitsio -lm
 
 map2fig_SOURCES = map2fig.c gopt.c