add phi0 parameter to sharp_make_gauss_geom_info()

This commit is contained in:
Martin Reinecke 2012-11-22 23:02:42 +01:00
parent 329f08b7b4
commit 6bbd0f46f5
9 changed files with 20 additions and 19 deletions

View file

@ -648,7 +648,7 @@ static int sharp_oracle (sharp_jobtype type, int spin, int ntrans)
ptrdiff_t npix=(ptrdiff_t)nrings*ppring; ptrdiff_t npix=(ptrdiff_t)nrings*ppring;
sharp_geom_info *tinfo; sharp_geom_info *tinfo;
sharp_make_gauss_geom_info (nrings, ppring, 1, ppring, &tinfo); sharp_make_gauss_geom_info (nrings, ppring, 0., 1, ppring, &tinfo);
ptrdiff_t nalms = ((mmax+1)*(mmax+2))/2 + (mmax+1)*(lmax-mmax); ptrdiff_t nalms = ((mmax+1)*(mmax+2))/2 + (mmax+1)*(lmax-mmax);
int ncomp = ntrans*((spin==0) ? 1 : 2); int ncomp = ntrans*((spin==0) ? 1 : 2);

View file

@ -94,7 +94,7 @@ static void check_sign_scale(void)
int nrings=lmax+1; int nrings=lmax+1;
int ppring=2*lmax+2; int ppring=2*lmax+2;
ptrdiff_t npix=(ptrdiff_t)nrings*ppring; ptrdiff_t npix=(ptrdiff_t)nrings*ppring;
sharp_make_gauss_geom_info (nrings, ppring, 1, ppring, &tinfo); sharp_make_gauss_geom_info (nrings, ppring, 0., 1, ppring, &tinfo);
/* flip theta to emulate the "old" Gaussian grid geometry */ /* flip theta to emulate the "old" Gaussian grid geometry */
for (int i=0; i<tinfo->npairs; ++i) for (int i=0; i<tinfo->npairs; ++i)
@ -247,7 +247,7 @@ int main(void)
int nrings=lmax+1; int nrings=lmax+1;
int ppring=2*lmax+2; int ppring=2*lmax+2;
ptrdiff_t npix=(ptrdiff_t)nrings*ppring; ptrdiff_t npix=(ptrdiff_t)nrings*ppring;
sharp_make_gauss_geom_info (nrings, ppring, 1, ppring, &tinfo); sharp_make_gauss_geom_info (nrings, ppring, 0., 1, ppring, &tinfo);
for (int nv=1; nv<=6; ++nv) for (int nv=1; nv<=6; ++nv)
for (int ntrans=1; ntrans<=6; ++ntrans) for (int ntrans=1; ntrans<=6; ++ntrans)
{ {

View file

@ -50,7 +50,7 @@ static void bench_sht (int spin, int nv, sharp_jobtype type,
int ppring=1024; int ppring=1024;
ptrdiff_t npix=(ptrdiff_t)nrings*ppring; ptrdiff_t npix=(ptrdiff_t)nrings*ppring;
sharp_geom_info *tinfo; sharp_geom_info *tinfo;
sharp_make_gauss_geom_info (nrings, ppring, 1, ppring, &tinfo); sharp_make_gauss_geom_info (nrings, ppring, 0., 1, ppring, &tinfo);
ptrdiff_t nalms = ((mmax+1)*(mmax+2))/2 + (mmax+1)*(lmax-mmax); ptrdiff_t nalms = ((mmax+1)*(mmax+2))/2 + (mmax+1)*(lmax-mmax);
int ncomp = ntrans*((spin==0) ? 1 : 2); int ncomp = ntrans*((spin==0) ? 1 : 2);

View file

@ -140,7 +140,7 @@ int main(int argc, char **argv)
{ {
int nrings=geom2=lmax+1; int nrings=geom2=lmax+1;
int ppring=atoi(argv[3]); int ppring=atoi(argv[3]);
sharp_make_gauss_geom_info (nrings, ppring, 1, ppring, &tinfo); sharp_make_gauss_geom_info (nrings, ppring, 0., 1, ppring, &tinfo);
} }
else if (strcmp(argv[1],"ecp")==0) else if (strcmp(argv[1],"ecp")==0)
{ {

View file

@ -173,15 +173,15 @@ static void makeweights (int bw, double *weights)
} }
} }
void sharp_make_gauss_geom_info (int nrings, int nphi, int stride_lon, void sharp_make_gauss_geom_info (int nrings, int nphi, double phi0,
int stride_lat, sharp_geom_info **geom_info) int stride_lon, int stride_lat, sharp_geom_info **geom_info)
{ {
const double pi=3.141592653589793238462643383279502884197; const double pi=3.141592653589793238462643383279502884197;
double *theta=RALLOC(double,nrings); double *theta=RALLOC(double,nrings);
double *weight=RALLOC(double,nrings); double *weight=RALLOC(double,nrings);
int *nph=RALLOC(int,nrings); int *nph=RALLOC(int,nrings);
double *phi0=RALLOC(double,nrings); double *phi0_=RALLOC(double,nrings);
ptrdiff_t *ofs=RALLOC(ptrdiff_t,nrings); ptrdiff_t *ofs=RALLOC(ptrdiff_t,nrings);
int *stride_=RALLOC(int,nrings); int *stride_=RALLOC(int,nrings);
@ -190,19 +190,19 @@ void sharp_make_gauss_geom_info (int nrings, int nphi, int stride_lon,
{ {
theta[m] = acos(-theta[m]); theta[m] = acos(-theta[m]);
nph[m]=nphi; nph[m]=nphi;
phi0[m]=0; phi0_[m]=phi0;
ofs[m]=(ptrdiff_t)m*stride_lat; ofs[m]=(ptrdiff_t)m*stride_lat;
stride_[m]=stride_lon; stride_[m]=stride_lon;
weight[m]*=2*pi/nphi; weight[m]*=2*pi/nphi;
} }
sharp_make_geom_info (nrings, nph, ofs, stride_, phi0, theta, NULL, weight, sharp_make_geom_info (nrings, nph, ofs, stride_, phi0_, theta, NULL, weight,
geom_info); geom_info);
DEALLOC(theta); DEALLOC(theta);
DEALLOC(weight); DEALLOC(weight);
DEALLOC(nph); DEALLOC(nph);
DEALLOC(phi0); DEALLOC(phi0_);
DEALLOC(ofs); DEALLOC(ofs);
DEALLOC(stride_); DEALLOC(stride_);
} }

View file

@ -53,12 +53,13 @@ void sharp_make_weighted_healpix_geom_info (int nside, int stride,
/*! Creates a geometry information describing a Gaussian map with \a nrings /*! Creates a geometry information describing a Gaussian map with \a nrings
iso-latitude rings and \a nphi pixels per ring. The azimuth of the first iso-latitude rings and \a nphi pixels per ring. The azimuth of the first
pixel in each ring is 0. The index difference between two adjacent pixels pixel in each ring is \a phi0 (in radians). The index difference between
in an iso-latitude ring is \a stride_lon, the index difference between the two adjacent pixels in an iso-latitude ring is \a stride_lon, the index
two start pixels in consecutive iso-latitude rings is \a stride_lat. difference between the two start pixels in consecutive iso-latitude rings
is \a stride_lat.
\ingroup geominfogroup */ \ingroup geominfogroup */
void sharp_make_gauss_geom_info (int nrings, int nphi, int stride_lon, void sharp_make_gauss_geom_info (int nrings, int nphi, double phi0,
int stride_lat, sharp_geom_info **geom_info); int stride_lon, int stride_lat, sharp_geom_info **geom_info);
/*! Creates a geometry information describing an ECP map with \a nrings /*! Creates a geometry information describing an ECP map with \a nrings
iso-latitude rings and \a nphi pixels per ring. The azimuth of the first iso-latitude rings and \a nphi pixels per ring. The azimuth of the first

View file

@ -190,7 +190,7 @@ int main(int argc, char **argv)
ptrdiff_t npix=0; ptrdiff_t npix=0;
int nrings=lmax+1; int nrings=lmax+1;
int ppring=atoi(argv[2]); int ppring=atoi(argv[2]);
sharp_make_gauss_geom_info (nrings, ppring, 1, ppring, &tinfo); sharp_make_gauss_geom_info (nrings, ppring, 0., 1, ppring, &tinfo);
reduce_geom_info(tinfo); reduce_geom_info(tinfo);
npix=get_npix(tinfo); npix=get_npix(tinfo);

View file

@ -211,7 +211,7 @@ int main(int argc, char **argv)
ptrdiff_t npix=(ptrdiff_t)nrings*ppring; ptrdiff_t npix=(ptrdiff_t)nrings*ppring;
printf("\nTesting Gaussian grid (%d rings, %d pixels/ring, %ld pixels)\n", printf("\nTesting Gaussian grid (%d rings, %d pixels/ring, %ld pixels)\n",
nrings,ppring,(long)npix); nrings,ppring,(long)npix);
sharp_make_gauss_geom_info (nrings, ppring, 1, ppring, &tinfo); sharp_make_gauss_geom_info (nrings, ppring, 0., 1, ppring, &tinfo);
check_accuracy(tinfo,lmax,lmax,npix,spin,ntrans,niter); check_accuracy(tinfo,lmax,lmax,npix,spin,ntrans,niter);
sharp_destroy_geom_info(tinfo); sharp_destroy_geom_info(tinfo);
} }

View file

@ -308,7 +308,7 @@ int main(int argc, char **argv)
if (mytask==0) if (mytask==0)
printf("\nTesting Gaussian grid (%d rings, %d pixels/ring, %ld pixels)\n", printf("\nTesting Gaussian grid (%d rings, %d pixels/ring, %ld pixels)\n",
nrings,ppring,(long)npix); nrings,ppring,(long)npix);
sharp_make_gauss_geom_info (nrings, ppring, 1, ppring, &tinfo); sharp_make_gauss_geom_info (nrings, ppring, 0., 1, ppring, &tinfo);
reduce_geom_info(tinfo); reduce_geom_info(tinfo);
npix=get_npix(tinfo); npix=get_npix(tinfo);
check_accuracy(tinfo,lmax,lmax,npix,spin,ntrans,niter); check_accuracy(tinfo,lmax,lmax,npix,spin,ntrans,niter);