Skip to content

Commit

Permalink
Fix healpix 'auto' tiling
Browse files Browse the repository at this point in the history
'auto' nside_tile is capped to nside
Add check that nside_tile <= nside in C++
  • Loading branch information
earosenberg committed Nov 7, 2024
1 parent b665511 commit bb392eb
Show file tree
Hide file tree
Showing 2 changed files with 7 additions and 2 deletions.
4 changes: 2 additions & 2 deletions python/proj/wcs.py
Original file line number Diff line number Diff line change
Expand Up @@ -766,15 +766,15 @@ def compute_nside_tile(self, assembly, nActivePerThread=5, nThreads=None):
"""
if self.nside_tile == 'auto':
## Estimate fsky
nside_tile0 = 4 # ntile = 192, for estimating fsky
nside_tile0 = min(4, self.nside) # ntile = 192, for estimating fsky
self.nside_tile = nside_tile0
nActive = len(self.get_active_tiles(assembly)['active_tiles'])
fsky = nActive / (12 * nside_tile0**2)
if nThreads is None:
nThreads = so3g.useful_info()['omp_num_threads']
# nside_tile is smallest power of 2 satisfying nTile >= nActivePerThread * nthread / fsky
self.nside_tile = int(2**np.ceil(0.5 * np.log2(nActivePerThread * nThreads / (12 * fsky))))
# print('Setting nside_tile =', self.nside_tile)
self.nside_tile = min(self.nside_tile, self.nside)
return self.nside_tile

def get_active_tiles(self, assembly, assign=False):
Expand Down
5 changes: 5 additions & 0 deletions src/Projection.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -596,6 +596,11 @@ class Pixelizor_Healpix<Tiled> {
npix_per_tile = npix / ntiles;
check_nside(nside); // check validity
check_nside(nside_tile);
if (nside_tile > nside) {
std::ostringstream err;
err << "Invalid nside_tile " << nside_tile << " > nside " << nside;
throw ValueError_exception(err.str());
}
}
~Pixelizor_Healpix() {};

Expand Down

0 comments on commit bb392eb

Please sign in to comment.