sharp.h 7.25 KB
 Martin Reinecke committed Jan 04, 2020 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 ``````/* * This file is part of libsharp2. * * libsharp2 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 2 of the License, or * (at your option) any later version. * * libsharp2 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 libsharp2; if not, write to the Free Software * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA */ /* libsharp2 is being developed at the Max-Planck-Institut fuer Astrophysik */ /*! \file sharp.h * Portable interface for the spherical transform library. * * Copyright (C) 2012-2019 Max-Planck-Society * \author Martin Reinecke \author Dag Sverre Seljebotn */ #ifndef SHARP_SHARP_H #define SHARP_SHARP_H `````` Martin Reinecke committed Jan 05, 2020 31 32 ``````#include #include `````` Martin Reinecke committed Jan 09, 2020 33 ``````#include `````` Martin Reinecke committed Jan 04, 2020 34 35 36 `````` /*! \internal Type holding all required information about a map geometry. */ `````` Martin Reinecke committed Jan 09, 2020 37 ``````struct sharp_geom_info `````` Martin Reinecke committed Jan 04, 2020 38 `````` { `````` Martin Reinecke committed Jan 12, 2020 39 40 41 42 43 44 45 46 47 48 49 `````` struct Tring { double theta, phi0, weight, cth, sth; ptrdiff_t ofs; int nph, stride; }; struct Tpair { Tring r1,r2; }; std::vector pair; `````` Martin Reinecke committed Jan 05, 2020 50 `````` int nphmax; `````` Martin Reinecke committed Jan 09, 2020 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 ``````/*! Creates a geometry information from a set of ring descriptions. All arrays passed to this function must have \a nrings elements. \param nrings the number of rings in the map \param nph the number of pixels in each ring \param ofs the index of the first pixel in each ring in the map array \param stride the stride between consecutive pixels \param phi0 the azimuth (in radians) of the first pixel in each ring \param theta the colatitude (in radians) of each ring \param wgt the pixel weight to be used for the ring in map2alm and adjoint map2alm transforms. Pass nullptr to use 1.0 as weight for all rings. */ sharp_geom_info(int nrings, const int *nph, const ptrdiff_t *ofs, const int *stride, const double *phi0, const double *theta, const double *wgt); `````` Martin Reinecke committed Jan 09, 2020 66 67 `````` void clear_map(double *map) const; void clear_map(float *map) const; `````` Martin Reinecke committed Jan 09, 2020 68 `````` }; `````` Martin Reinecke committed Jan 04, 2020 69 70 71 72 73 74 `````` /*! \defgroup almgroup Helpers for dealing with a_lm */ /*! \{ */ /*! \internal Helper type for index calculation in a_lm arrays. */ `````` Martin Reinecke committed Jan 09, 2020 75 ``````struct sharp_alm_info `````` Martin Reinecke committed Jan 04, 2020 76 77 78 79 80 81 `````` { /*! Maximum \a l index of the array */ int lmax; /*! Number of different \a m values in this object */ int nm; /*! Array with \a nm entries containing the individual m values */ `````` Martin Reinecke committed Jan 05, 2020 82 `````` std::vector mval; `````` Martin Reinecke committed Jan 04, 2020 83 84 `````` /*! Array with \a nm entries containing the (hypothetical) indices of the coefficients with quantum numbers 0,\a mval[i] */ `````` Martin Reinecke committed Jan 05, 2020 85 `````` std::vector mvstart; `````` Martin Reinecke committed Jan 04, 2020 86 87 88 89 90 91 92 93 94 95 96 `````` /*! Stride between a_lm and a_(l+1),m */ ptrdiff_t stride; /*! Creates an a_lm data structure from the following parameters: \param lmax maximum \a l quantum number (>=0) \param mmax maximum \a m quantum number (0<= \a mmax <= \a lmax) \param stride the stride between entries with identical \a m, and \a l differing by 1. \param mstart the index of the (hypothetical) coefficient with the quantum numbers 0,\a m. Must have \a mmax+1 entries. */ `````` Martin Reinecke committed Jan 09, 2020 97 98 `````` sharp_alm_info(int lmax, int mmax, int stride, const ptrdiff_t *mstart); `````` Martin Reinecke committed Jan 04, 2020 99 100 101 102 103 104 105 106 107 ``````/*! Creates an a_lm data structure which from the following parameters: \param lmax maximum \a l quantum number (\a >=0) \param nm number of different \a m (\a 0<=nm<=lmax+1) \param stride the stride between entries with identical \a m, and \a l differing by 1. \param mval array with \a nm entries containing the individual m values \param mvstart array with \a nm entries containing the (hypothetical) indices of the coefficients with the quantum numbers 0,\a mval[i] */ `````` Martin Reinecke committed Jan 09, 2020 108 `````` sharp_alm_info (int lmax, int nm, int stride, const int *mval, `````` Martin Reinecke committed Jan 09, 2020 109 `````` const ptrdiff_t *mvstart); `````` Martin Reinecke committed Jan 04, 2020 110 111 112 113 ``````/*! Returns the index of the coefficient with quantum numbers \a l, \a mval[mi]. \note for a \a sharp_alm_info generated by sharp_make_alm_info() this is the index for the coefficient with the quantum numbers \a l, \a mi. */ `````` Martin Reinecke committed Jan 09, 2020 114 115 116 `````` ptrdiff_t index (int l, int mi); }; `````` Martin Reinecke committed Jan 04, 2020 117 118 119 120 121 122 123 124 125 126 127 ``````/*! \} */ /*! \defgroup geominfogroup Functions for dealing with geometry information */ /*! \{ */ /*! \} */ /*! \defgroup lowlevelgroup Low-level libsharp2 SHT interface */ /*! \{ */ /*! Enumeration of SHARP job types. */ `````` Martin Reinecke committed Jan 09, 2020 128 ``````enum sharp_jobtype { SHARP_YtW=0, /*!< analysis */ `````` Martin Reinecke committed Jan 04, 2020 129 130 131 132 133 134 `````` SHARP_MAP2ALM=SHARP_YtW, /*!< analysis */ SHARP_Y=1, /*!< synthesis */ SHARP_ALM2MAP=SHARP_Y, /*!< synthesis */ SHARP_Yt=2, /*!< adjoint synthesis */ SHARP_WY=3, /*!< adjoint analysis */ SHARP_ALM2MAP_DERIV1=4 /*!< synthesis of first derivatives */ `````` Martin Reinecke committed Jan 09, 2020 135 `````` }; `````` Martin Reinecke committed Jan 04, 2020 136 137 `````` /*! Job flags */ `````` Martin Reinecke committed Jan 09, 2020 138 ``````enum sharp_jobflags { SHARP_DP = 1<<4, `````` Martin Reinecke committed Jan 04, 2020 139 140 141 142 143 `````` /*!< map and a_lm are in double precision */ SHARP_ADD = 1<<5, /*!< results are added to the output arrays, instead of overwriting them */ SHARP_USE_WEIGHTS = 1<<20, /* internal use only */ `````` Martin Reinecke committed Jan 09, 2020 144 `````` }; `````` Martin Reinecke committed Jan 04, 2020 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 `````` /*! Performs a libsharp2 SHT job. The interface deliberately does not use the C99 "complex" data type, in order to be callable from C89 and C++. \param type the type of SHT \param spin the spin of the quantities to be transformed \param alm contains pointers to the a_lm coefficients. If \a spin==0, alm[0] points to the a_lm of the SHT. If \a spin>0, alm[0] and alm[1] point to the two a_lm sets of the SHT. The exact data type of \a alm depends on whether the SHARP_DP flag is set. \param map contains pointers to the maps. If \a spin==0, map[0] points to the map of the SHT. If \a spin>0, or \a type is SHARP_ALM2MAP_DERIV1, map[0] and map[1] point to the two maps of the SHT. The exact data type of \a map depends on whether the SHARP_DP flag is set. \param geom_info A \c sharp_geom_info object compatible with the provided \a map arrays. \param alm_info A \c sharp_alm_info object compatible with the provided \a alm arrays. All \c m values from 0 to some \c mmax<=lmax must be present exactly once. \param flags See sharp_jobflags. In particular, if SHARP_DP is set, then \a alm is expected to have the type "complex double **" and \a map is expected to have the type "double **"; otherwise, the expected types are "complex float **" and "float **", respectively. `````` Martin Reinecke committed Jan 09, 2020 167 `````` \param time If not nullptr, the wall clock time required for this SHT `````` Martin Reinecke committed Jan 04, 2020 168 `````` (in seconds) will be written here. `````` Martin Reinecke committed Jan 09, 2020 169 `````` \param opcnt If not nullptr, a conservative estimate of the total floating point `````` Martin Reinecke committed Jan 04, 2020 170 171 `````` operation count for this SHT will be written here. */ void sharp_execute (sharp_jobtype type, int spin, void *alm, void *map, `````` Martin Reinecke committed Jan 09, 2020 172 `````` const sharp_geom_info &geom_info, const sharp_alm_info &alm_info, `````` Martin Reinecke committed Jan 04, 2020 173 174 175 176 177 178 179 180 181 182 183 184 `````` int flags, double *time, unsigned long long *opcnt); void sharp_set_chunksize_min(int new_chunksize_min); void sharp_set_nchunks_max(int new_nchunks_max); /*! \} */ int sharp_get_mlim (int lmax, int spin, double sth, double cth); int sharp_veclen(void); const char *sharp_architecture(void); #endif``````