LevelS SHT library  3.50
sharp_lowlevel.h
Go to the documentation of this file.
1 /*
2  * This file is part of libsharp.
3  *
4  * libsharp is free software; you can redistribute it and/or modify
5  * it under the terms of the GNU General Public License as published by
6  * the Free Software Foundation; either version 2 of the License, or
7  * (at your option) any later version.
8  *
9  * libsharp is distributed in the hope that it will be useful,
10  * but WITHOUT ANY WARRANTY; without even the implied warranty of
11  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12  * GNU General Public License for more details.
13  *
14  * You should have received a copy of the GNU General Public License
15  * along with libsharp; if not, write to the Free Software
16  * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
17  */
18 
19 /*
20  * libsharp is being developed at the Max-Planck-Institut fuer Astrophysik
21  * and financially supported by the Deutsches Zentrum fuer Luft- und Raumfahrt
22  * (DLR).
23  */
24 
25 /*! \file sharp_lowlevel.h
26  * Low-level, portable interface for the spherical transform library.
27  *
28  * Copyright (C) 2012-2013 Max-Planck-Society
29  * \author Martin Reinecke \author Dag Sverre Seljebotn
30  */
31 
32 #ifndef PLANCK_SHARP_LOWLEVEL_H
33 #define PLANCK_SHARP_LOWLEVEL_H
34 
35 #include <stddef.h>
36 
37 #ifdef __cplusplus
38 extern "C" {
39 #endif
40 
41 /*! \internal
42  Helper type containing information about a single ring. */
43 typedef struct
44  {
45  double theta, phi0, weight, cth, sth;
46  ptrdiff_t ofs;
47  int nph, stride;
48  } sharp_ringinfo;
49 
50 /*! \internal
51  Helper type containing information about a pair of rings with colatitudes
52  symmetric around the equator. */
53 typedef struct
54  {
55  sharp_ringinfo r1,r2;
56  } sharp_ringpair;
57 
58 /*! \internal
59  Type holding all required information about a map geometry. */
60 typedef struct
61  {
62  sharp_ringpair *pair;
63  int npairs, nphmax;
64  } sharp_geom_info;
65 
66 /*! \defgroup almgroup Helpers for dealing with a_lm */
67 /*! \{ */
68 
69 /*! \internal
70  Helper type for index calculation in a_lm arrays. */
71 typedef struct
72  {
73  /*! Maximum \a l index of the array */
74  int lmax;
75  /*! Number of different \a m values in this object */
76  int nm;
77  /*! Array with \a nm entries containing the individual m values */
78  int *mval;
79  /*! Combination of flags from sharp_almflags */
80  int flags;
81  /*! Array with \a nm entries containing the (hypothetical) indices of
82  the coefficients with quantum numbers 0,\a mval[i] */
83  ptrdiff_t *mvstart;
84  /*! Stride between a_lm and a_(l+1),m */
85  ptrdiff_t stride;
86  } sharp_alm_info;
87 
88 /*! alm_info flags */
89 typedef enum { SHARP_PACKED = 1,
90  /*!< m=0-coefficients are packed so that the (zero) imaginary part is
91  not present. mvstart is in units of *real* float/double for all
92  m; stride is in units of reals for m=0 and complex for m!=0 */
94  /*!< Use the real spherical harmonic convention. For
95  m==0, the alm are treated exactly the same as in
96  the complex case. For m!=0, alm[i] represent a
97  pair (+abs(m), -abs(m)) instead of (real, imag),
98  and the coefficients are scaled by a factor of
99  sqrt(2) relative to the complex case. In other
100  words, (sqrt(.5) * alm[i]) recovers the
101  corresponding complex coefficient (when accessed
102  as complex).
103  */
104  } sharp_almflags;
105 
106 
107 
108 /*! Creates an a_lm data structure from the following parameters:
109  \param lmax maximum \a l quantum number (>=0)
110  \param mmax maximum \a m quantum number (0<= \a mmax <= \a lmax)
111  \param stride the stride between entries with identical \a m, and \a l
112  differing by 1.
113  \param mstart the index of the (hypothetical) coefficient with the
114  quantum numbers 0,\a m. Must have \a mmax+1 entries.
115  \param alm_info will hold a pointer to the newly created data structure
116  */
117 void sharp_make_alm_info (int lmax, int mmax, int stride,
118  const ptrdiff_t *mstart, sharp_alm_info **alm_info);
119 /*! Creates an a_lm data structure which from the following parameters:
120  \param lmax maximum \a l quantum number (\a >=0)
121  \param nm number of different \a m (\a 0<=nm<=lmax+1)
122  \param stride the stride between entries with identical \a m, and \a l
123  differing by 1.
124  \param mval array with \a nm entries containing the individual m values
125  \param mvstart array with \a nm entries containing the (hypothetical)
126  indices of the coefficients with the quantum numbers 0,\a mval[i]
127  \param flags a combination of sharp_almflags (pass 0 unless you know you need this)
128  \param alm_info will hold a pointer to the newly created data structure
129  */
130 void sharp_make_general_alm_info (int lmax, int nm, int stride, const int *mval,
131  const ptrdiff_t *mvstart, int flags, sharp_alm_info **alm_info);
132 /*! Returns the index of the coefficient with quantum numbers \a l,
133  \a mval[mi].
134  \note for a \a sharp_alm_info generated by sharp_make_alm_info() this is
135  the index for the coefficient with the quantum numbers \a l, \a mi. */
136 ptrdiff_t sharp_alm_index (const sharp_alm_info *self, int l, int mi);
137 /*! Returns the number of alm coefficients described by \a self. If the SHARP_PACKED
138  flag is set, this is number of "real" coeffecients (for m < 0 and m >= 0),
139  otherwise it is the number of complex coefficients (with m>=0). */
140 ptrdiff_t sharp_alm_count(const sharp_alm_info *self);
141 /*! Deallocates the a_lm info object. */
142 void sharp_destroy_alm_info (sharp_alm_info *info);
143 
144 /*! \} */
145 
146 /*! \defgroup geominfogroup Functions for dealing with geometry information */
147 /*! \{ */
148 
149 /*! Creates a geometry information from a set of ring descriptions.
150  All arrays passed to this function must have \a nrings elements.
151  \param nrings the number of rings in the map
152  \param nph the number of pixels in each ring
153  \param ofs the index of the first pixel in each ring in the map array
154  \param stride the stride between consecutive pixels
155  \param phi0 the azimuth (in radians) of the first pixel in each ring
156  \param theta the colatitude (in radians) of each ring
157  \param wgt the pixel weight to be used for the ring in map2alm
158  and adjoint map2alm transforms.
159  Pass NULL to use 1.0 as weight for all rings.
160  \param geom_info will hold a pointer to the newly created data structure
161  */
162 void sharp_make_geom_info (int nrings, const int *nph, const ptrdiff_t *ofs,
163  const int *stride, const double *phi0, const double *theta,
164  const double *wgt, sharp_geom_info **geom_info);
165 
166 /*! Counts the number of grid points needed for (the local part of) a map described
167  by \a info.
168  */
169 ptrdiff_t sharp_map_size(const sharp_geom_info *info);
170 
171 /*! Deallocates the geometry information in \a info. */
172 void sharp_destroy_geom_info (sharp_geom_info *info);
173 
174 /*! \} */
175 
176 /*! \defgroup lowlevelgroup Low-level libsharp SHT interface */
177 /*! \{ */
178 
179 /*! Enumeration of SHARP job types. */
180 typedef enum { SHARP_YtW=0, /*!< analysis */
181  SHARP_MAP2ALM=SHARP_YtW, /*!< analysis */
182  SHARP_Y=1, /*!< synthesis */
183  SHARP_ALM2MAP=SHARP_Y, /*!< synthesis */
184  SHARP_Yt=2, /*!< adjoint synthesis */
185  SHARP_WY=3, /*!< adjoint analysis */
186  SHARP_ALM2MAP_DERIV1=4 /*!< synthesis of first derivatives */
187  } sharp_jobtype;
188 
189 /*! Job flags */
190 typedef enum { SHARP_DP = 1<<4,
191  /*!< map and a_lm are in double precision */
192  SHARP_ADD = 1<<5,
193  /*!< results are added to the output arrays, instead of
194  overwriting them */
195 
196  /* NOTE: SHARP_REAL_HARMONICS, 1<<6, is also available in sharp_jobflags,
197  but its use here is deprecated in favor of having it in the sharp_alm_info */
198 
199  SHARP_NO_FFT = 1<<7,
200 
201  SHARP_USE_WEIGHTS = 1<<20, /* internal use only */
202  SHARP_NO_OPENMP = 1<<21, /* internal use only */
203  SHARP_NVMAX = (1<<4)-1 /* internal use only */
204  } sharp_jobflags;
205 
206 /*! Performs a libsharp SHT job. The interface deliberately does not use
207  the C99 "complex" data type, in order to be callable from C89 and C++.
208  \param type the type of SHT
209  \param spin the spin of the quantities to be transformed
210  \param alm contains pointers to the a_lm coefficients. If \a spin==0,
211  alm[0] points to the a_lm of the first SHT, alm[1] to those of the second
212  etc. If \a spin>0, alm[0] and alm[1] point to the a_lm of the first SHT,
213  alm[2] and alm[3] to those of the second, etc. The exact data type of \a alm
214  depends on whether the SHARP_DP flag is set.
215  \param map contains pointers to the maps. If \a spin==0,
216  map[0] points to the map of the first SHT, map[1] to that of the second
217  etc. If \a spin>0, or \a type is SHARP_ALM2MAP_DERIV1, map[0] and map[1]
218  point to the maps of the first SHT, map[2] and map[3] to those of the
219  second, etc. The exact data type of \a map depends on whether the SHARP_DP
220  flag is set.
221  \param geom_info A \c sharp_geom_info object compatible with the provided
222  \a map arrays.
223  \param alm_info A \c sharp_alm_info object compatible with the provided
224  \a alm arrays. All \c m values from 0 to some \c mmax<=lmax must be present
225  exactly once.
226  \param ntrans the number of simultaneous SHTs
227  \param flags See sharp_jobflags. In particular, if SHARP_DP is set, then
228  \a alm is expected to have the type "complex double **" and \a map is
229  expected to have the type "double **"; otherwise, the expected
230  types are "complex float **" and "float **", respectively.
231  \param time If not NULL, the wall clock time required for this SHT
232  (in seconds) will be written here.
233  \param opcnt If not NULL, a conservative estimate of the total floating point
234  operation count for this SHT will be written here. */
235 void sharp_execute (sharp_jobtype type, int spin, void *alm, void *map,
236  const sharp_geom_info *geom_info, const sharp_alm_info *alm_info, int ntrans,
237  int flags, double *time, unsigned long long *opcnt);
238 
239 void sharp_set_chunksize_min(int new_chunksize_min);
240 void sharp_set_nchunks_max(int new_nchunks_max);
241 
242 /*! \} */
243 
244 #ifdef __cplusplus
245 }
246 #endif
247 
248 #endif
void sharp_execute(sharp_jobtype type, int spin, void *alm, void *map, const sharp_geom_info *geom_info, const sharp_alm_info *alm_info, int ntrans, int flags, double *time, unsigned long long *opcnt)
Definition: sharp.c:950
void sharp_destroy_geom_info(sharp_geom_info *info)
Definition: sharp.c:259
sharp_jobtype
void sharp_make_general_alm_info(int lmax, int nm, int stride, const int *mval, const ptrdiff_t *mvstart, int flags, sharp_alm_info **alm_info)
Definition: sharp.c:143
void sharp_destroy_alm_info(sharp_alm_info *info)
Definition: sharp.c:191
sharp_jobflags
ptrdiff_t sharp_map_size(const sharp_geom_info *info)
Definition: sharp.c:248
void sharp_make_geom_info(int nrings, const int *nph, const ptrdiff_t *ofs, const int *stride, const double *phi0, const double *theta, const double *wgt, sharp_geom_info **geom_info)
Definition: sharp.c:198
sharp_almflags
ptrdiff_t sharp_alm_index(const sharp_alm_info *self, int l, int mi)
Definition: sharp.c:171
void sharp_make_alm_info(int lmax, int mmax, int stride, const ptrdiff_t *mstart, sharp_alm_info **alm_info)
Definition: sharp.c:161
ptrdiff_t sharp_alm_count(const sharp_alm_info *self)
Definition: sharp.c:178

Generated on Mon Dec 10 2018 10:24:20 for LevelS SHT library