shtns.h

Go to the documentation of this file.

/*
 * Copyright (c) 2010-2021 Centre National de la Recherche Scientifique.
 * written by Nathanael Schaeffer (CNRS, ISTerre, Grenoble, France).
 * 
 * nathanael.schaeffer@univ-grenoble-alpes.fr
 * 
 * This software is governed by the CeCILL license under French law and
 * abiding by the rules of distribution of free software. You can use,
 * modify and/or redistribute the software under the terms of the CeCILL
 * license as circulated by CEA, CNRS and INRIA at the following URL
 * "http://www.cecill.info".
 * 
 * The fact that you are presently reading this means that you have had
 * knowledge of the CeCILL license and that you accept its terms.
 * 
 */
 
#ifdef __cplusplus
        #include <complex>
        typedef std::complex<double> cplx;              
        typedef std::complex<float> cplx_f;             
#else
        #include <complex.h>
        typedef _Complex double cplx;                   
        typedef _Complex float cplx_f;                  
#endif
 
#ifdef __cplusplus
extern "C" {
#endif /* __cplusplus */
 
#define SHTNS_INTERFACE 0x30602
 
typedef struct shtns_info* shtns_cfg;
 
typedef struct shtns_rot_* shtns_rot;
 
enum shtns_norm {
        sht_orthonormal,        
        sht_fourpi,                     
        sht_schmidt,            
        sht_for_rotations       
};
#define SHT_NO_CS_PHASE (256*4)         
#define SHT_REAL_NORM (256*8)           
 
enum shtns_type {
        sht_gauss,              
        sht_auto,               
        sht_reg_fast,   
        sht_reg_dct,    
        sht_quick_init, 
        sht_reg_poles,  
        sht_gauss_fly   
};
#define SHT_NATIVE_LAYOUT 0                     
#define SHT_THETA_CONTIGUOUS 256        
#define SHT_PHI_CONTIGUOUS (256*2)      
#define SHT_SOUTH_POLE_FIRST (256*32)   
 
#define SHT_SCALAR_ONLY (256*16)        
#define SHT_LOAD_SAVE_CFG (256*64)      
#define SHT_ALLOW_GPU (256*128)         
#define SHT_ALLOW_PADDING (256*256)     
#define SHT_ROBERT_FORM (256*512)       
#define SHT_FP32 (256*1024)                     
 
 
#ifndef SHTNS_PRIVATE
struct shtns_info {             // allow read-only access to some data (useful for optimization and helper macros)
        const unsigned int nlm;                 
        const unsigned short lmax;              
        const unsigned short mmax;              
        const unsigned short mres;              
        const unsigned short nlat_2;    
        const unsigned int nlat;                
        const unsigned int nphi;                
        const unsigned int nspat;               
        const unsigned short *const li; 
        const unsigned short *const mi; 
        const double *const ct;                 
        const double *const st;                 
        const unsigned int nlat_padded; 
        const unsigned int nlm_cplx;    
};
#endif
 
// MACROS //
 
 
#define LiM(shtns, l,im) ( (((im)*(2*shtns->lmax + 2 - ((im)+1)*shtns->mres))>>1) + (l) )
#define LM(shtns, l,m) ( (((((unsigned short)(m))/shtns->mres)*(2*shtns->lmax + 2 - ((m)+shtns->mres)))>>1) + (l) )
#define LM_LOOP( shtns, action ) { int lm=0; do { action } while(++lm < shtns->nlm); }
#define LM_L_LOOP( shtns, action ) { int lm=0; do { int l=shtns->li[lm]; action } while(++lm < shtns->nlm); }
#define LM_cplx(shtns, l, m) ( ((l) <= shtns->mmax) ? (l)*((l)+1)+(m) : shtns->mmax*(2*(l) - shtns->mmax) + (l)+(m) )
 
 
#define NSPAT_ALLOC(shtns) (shtns->nspat)
 
// HELPER MACROS //
 
#define PHI_DEG(shtns, ip) (360./(shtns->nphi*shtns->mres))*(ip)
#define PHI_RAD(shtns, ip) (2.*M_PI/(shtns->nphi*shtns->mres))*(ip)
 
 
// FUNCTIONS //
 
long nlm_calc(long lmax, long mmax, long mres);
 
long nlm_cplx_calc(long lmax, long mmax, long mres);
 
void shtns_verbose(int);                        
void shtns_print_version(void);         
const char* shtns_get_build_info(void); 
 
#ifndef SWIG
 
void shtns_print_cfg(shtns_cfg);        
 
 
shtns_cfg shtns_init(enum shtns_type flags, int lmax, int mmax, int mres, int nlat, int nphi);
shtns_cfg shtns_create(int lmax, int mmax, int mres, enum shtns_norm norm);
int shtns_set_grid(shtns_cfg, enum shtns_type flags, double eps, int nlat, int nphi);
int shtns_set_grid_auto(shtns_cfg, enum shtns_type flags, double eps, int nl_order, int *nlat, int *nphi);
shtns_cfg shtns_create_with_grid(shtns_cfg, int mmax, int nofft);
int shtns_use_threads(int num_threads);
int shtns_set_batch(shtns_cfg shtns, int howmany, long spec_dist);
 
void shtns_reset(void);                         
void shtns_destroy(shtns_cfg);          
void shtns_unset_grid(shtns_cfg);       
 
void shtns_robert_form(shtns_cfg, int robert);
 
void* shtns_malloc(size_t bytes);       
void shtns_free(void* p);                       
 
 
 
double sh00_1(shtns_cfg);       
double sh10_ct(shtns_cfg);      
double sh11_st(shtns_cfg);      
double shlm_e1(shtns_cfg, int l, int m);                
int shtns_gauss_wts(shtns_cfg, double *wts);
 
 
 
void SH_Zrotate(shtns_cfg, cplx *Qlm, double alpha, cplx *Rlm);
void SH_Yrotate(shtns_cfg, cplx *Qlm, double alpha, cplx *Rlm);
void SH_Yrotate90(shtns_cfg, cplx *Qlm, cplx *Rlm);
void SH_Xrotate90(shtns_cfg, cplx *Qlm, cplx *Rlm);
 
 
shtns_rot shtns_rotation_create(const int lmax, const int mmax, int norm);
void shtns_rotation_destroy(shtns_rot r);
void shtns_rotation_set_angles_ZYZ(shtns_rot r, double alpha, double beta, double gamma);
void shtns_rotation_set_angles_ZXZ(shtns_rot r, double alpha, double beta, double gamma);
void shtns_rotation_set_angle_axis(shtns_rot r, double theta, double Vx, double Vy, double Vz);
int shtns_rotation_wigner_d_matrix(shtns_rot r, const int l, double* mx);
void shtns_rotation_apply_cplx(shtns_rot r, cplx* Zlm, cplx* Rlm);
void shtns_rotation_apply_real(shtns_rot r, cplx* Qlm, cplx* Rlm);
 
 
 
int legendre_sphPlm_array(shtns_cfg shtns, const int lmax, const int im, const double x, double *yl);
int legendre_sphPlm_deriv_array(shtns_cfg shtns, const int lmax, const int im, const double x, const double sint, double *yl, double *dyl);
 
 
void mul_ct_matrix(shtns_cfg, double* mx);
void st_dt_matrix(shtns_cfg, double* mx);
void SH_mul_mx(shtns_cfg, double* mx, cplx *Qlm, cplx *Rlm);
 
 
 
void spat_to_SH(shtns_cfg shtns, double *Vr, cplx *Qlm);
void SH_to_spat(shtns_cfg shtns, cplx *Qlm, double *Vr);
void SH_to_spat_cplx(shtns_cfg shtns, cplx *alm, cplx *z);
void spat_cplx_to_SH(shtns_cfg shtns, cplx *z, cplx *alm);
 
 
void spat_to_SHsphtor(shtns_cfg, double *Vt, double *Vp, cplx *Slm, cplx *Tlm);
void SHsphtor_to_spat(shtns_cfg, cplx *Slm, cplx *Tlm, double *Vt, double *Vp);
void SHsph_to_spat(shtns_cfg, cplx *Slm, double *Vt, double *Vp);
void SHtor_to_spat(shtns_cfg, cplx *Tlm, double *Vt, double *Vp);
 
void spat_cplx_to_SHsphtor(shtns_cfg, cplx *Vt, cplx *Vp, cplx *Slm, cplx *Tlm);
void SHsphtor_to_spat_cplx(shtns_cfg, cplx *Slm, cplx *Tlm, cplx *Vt, cplx *Vp);
#define SH_to_grad_spat(shtns, S,Gt,Gp) SHsph_to_spat(shtns, S, Gt, Gp)
 
 
void spat_to_SHqst(shtns_cfg, double *Vr, double *Vt, double *Vp, cplx *Qlm, cplx *Slm, cplx *Tlm);
void SHqst_to_spat(shtns_cfg, cplx *Qlm, cplx *Slm, cplx *Tlm, double *Vr, double *Vt, double *Vp);
 
void spat_cplx_to_SHqst(shtns_cfg, cplx *Vr, cplx *Vt, cplx *Vp, cplx *Qlm, cplx *Slm, cplx *Tlm);
void SHqst_to_spat_cplx(shtns_cfg, cplx *Qlm, cplx *Slm, cplx *Tlm, cplx *Vr, cplx *Vt, cplx *Vp);
 
void spat_to_SH_l(shtns_cfg, double *Vr, cplx *Qlm, int ltr);
void SH_to_spat_l(shtns_cfg, cplx *Qlm, double *Vr, int ltr);
 
void SHsphtor_to_spat_l(shtns_cfg, cplx *Slm, cplx *Tlm, double *Vt, double *Vp, int ltr);
void SHsph_to_spat_l(shtns_cfg, cplx *Slm, double *Vt, double *Vp, int ltr);
void SHtor_to_spat_l(shtns_cfg, cplx *Tlm, double *Vt, double *Vp, int ltr);
void spat_to_SHsphtor_l(shtns_cfg, double *Vt, double *Vp, cplx *Slm, cplx *Tlm, int ltr);
 
void spat_to_SHqst_l(shtns_cfg, double *Vr, double *Vt, double *Vp, cplx *Qlm, cplx *Slm, cplx *Tlm, int ltr);
void SHqst_to_spat_l(shtns_cfg, cplx *Qlm, cplx *Slm, cplx *Tlm, double *Vr, double *Vt, double *Vp, int ltr);
#define SH_to_grad_spat_l(shtns, S,Gt,Gp,ltr) SHsph_to_spat_l(shtns, S, Gt, Gp, ltr)
 
void spat_to_SH_ml(shtns_cfg, int im, cplx *Vr, cplx *Ql, int ltr);
void SH_to_spat_ml(shtns_cfg, int im, cplx *Ql, cplx *Vr, int ltr);
 
void spat_to_SHsphtor_ml(shtns_cfg, int im, cplx *Vt, cplx *Vp, cplx *Sl, cplx *Tl, int ltr);
void SHsphtor_to_spat_ml(shtns_cfg, int im, cplx *Sl, cplx *Tl, cplx *Vt, cplx *Vp, int ltr);
void SHsph_to_spat_ml(shtns_cfg, int im, cplx *Sl, cplx *Vt, cplx *Vp, int ltr);
void SHtor_to_spat_ml(shtns_cfg, int im, cplx *Tl, cplx *Vt, cplx *Vp, int ltr);
 
void spat_to_SHqst_ml(shtns_cfg, int im, cplx *Vr, cplx *Vt, cplx *Vp, cplx *Ql, cplx *Sl, cplx *Tl, int ltr);
void SHqst_to_spat_ml(shtns_cfg, int im, cplx *Ql, cplx *Sl, cplx *Tl, cplx *Vr, cplx *Vt, cplx *Vp, int ltr);
#define SH_to_grad_spat_ml(shtns, im, S,Gt,Gp,ltr) SHsph_to_spat_ml(shtns, im, S, Gt, Gp, ltr)
 
 
double SH_to_point(shtns_cfg, cplx *Qlm, double cost, double phi);
void SH_to_point_cplx(shtns_cfg, cplx *alm, double cost, double phi, cplx* z_out);      
void SH_to_grad_point(shtns_cfg, cplx *DrSlm, cplx *Slm,
                                        double cost, double phi, double *vr, double *vt, double *vp);
void SHqst_to_point(shtns_cfg, cplx *Qlm, cplx *Slm, cplx *Tlm,
                                        double cost, double phi, double *vr, double *vt, double *vp);
 
void SH_to_lat(shtns_cfg shtns, cplx *Qlm, double cost,
                                        double *vr, int nphi, int ltr, int mtr);
void SHqst_to_lat(shtns_cfg, cplx *Qlm, cplx *Slm, cplx *Tlm, double cost,
                                        double *vr, double *vt, double *vp, int nphi, int ltr, int mtr);
 
 
#endif
 
#ifdef __cplusplus
}
#endif /* __cplusplus */