add MW geometry

libsharp/sharp_geomhelpers.c

@@ -291,3 +291,34 @@ void sharp_make_hw_geom_info (int nrings, int ppring, double phi0,
   DEALLOC(ofs);
   DEALLOC(stride_);
   }
+
+void sharp_make_mw_geom_info (int nrings, int ppring, double phi0,
+  int stride_lon, int stride_lat, sharp_geom_info **geom_info)
+  {
+  const double pi=3.141592653589793238462643383279502884197;
+
+  double *theta=RALLOC(double,nrings);
+  int *nph=RALLOC(int,nrings);
+  double *phi0_=RALLOC(double,nrings);
+  ptrdiff_t *ofs=RALLOC(ptrdiff_t,nrings);
+  int *stride_=RALLOC(int,nrings);
+
+  for (int m=0; m<nrings; ++m)
+    {
+    theta[m]=pi*(2.*m+1.)/(2.*nrings-1.);
+    if (theta[m]>pi-1e-15) theta[m]=pi-1e-15;
+    nph[m]=ppring;
+    phi0_[m]=phi0;
+    ofs[m]=(ptrdiff_t)m*stride_lat;
+    stride_[m]=stride_lon;
+    }
+
+  sharp_make_geom_info (nrings, nph, ofs, stride_, phi0_, theta, NULL, NULL,
+    geom_info);
+
+  DEALLOC(theta);
+  DEALLOC(nph);
+  DEALLOC(phi0_);
+  DEALLOC(ofs);
+  DEALLOC(stride_);
+  }

libsharp/sharp_geomhelpers.h

@@ -87,11 +87,27 @@ void sharp_make_ecp_geom_info (int nrings, int nphi, double phi0,
       longitude.
     \note The sphere is pixelized in a way that the colatitude of the first ring
       is \a 0 and that of the last ring is \a pi.
-    \note This is the grid used by Huffenberger & Wandelt 2010
+    \note This is the grid used by Huffenberger & Wandelt 2010.
     \ingroup geominfogroup */
 void sharp_make_hw_geom_info (int nrings, int ppring, double phi0,
   int stride_lon, int stride_lat, sharp_geom_info **geom_info);
 
+/*! Creates a geometry information describing a map with \a nrings
+    iso-latitude rings and \a nphi pixels per ring. The azimuth of the first
+    pixel in each ring is \a phi0 (in radians). The index difference between
+    two adjacent pixels in an iso-latitude ring is \a stride_lon, the index
+    difference between the two start pixels in consecutive iso-latitude rings
+    is \a stride_lat.
+    \note The spacing of pixel centers is equidistant in colatitude and
+      longitude.
+    \note The sphere is pixelized in a way that the colatitude of the first ring
+      is \a pi/(2*nrings-1) and that of the last ring is \a pi.
+    \note This is the grid introduced by McEwen & Wiaux 2011.
+    \note This function does \e not define any quadrature weights.
+    \ingroup geominfogroup */
+void sharp_make_mw_geom_info (int nrings, int ppring, double phi0,
+  int stride_lon, int stride_lat, sharp_geom_info **geom_info);
+
 #ifdef __cplusplus
 }
 #endif