Skip to content

Commit

Permalink
smurf: Only switch off wcs wrapping for tiles that straddle RA=12h
Browse files Browse the repository at this point in the history
  • Loading branch information
David Berry committed Sep 30, 2014
1 parent 26b29d3 commit 56bf8fa
Showing 1 changed file with 24 additions and 6 deletions.
30 changes: 24 additions & 6 deletions applications/smurf/libsmf/smf_jsadicer.c
Original file line number Diff line number Diff line change
Expand Up @@ -101,6 +101,10 @@
* the bounding box of each tile in the output NDF. This is to
* avoid tailes that span RA=12h being wrapped round to the far
* side of the projection, thus producing a huge box.
* 30-SEP-2014 (DSB):
* Only disable lon/lat wrapping for tiles that span RA=12h. The
* other tiles in the split HPX facet require lon/lat wrapping to
* be on, in order to work properly.
* {enter_further_changes_here}
* Copyright:
Expand Down Expand Up @@ -481,25 +485,39 @@ void smf_jsadicer( int indf, const char *base, int trim, smf_inst_t instrument,

/* Next job is to find the pixel bounds of the output NDF to create
which will hold data for the current tile. First map the pixel bounds
of the whole tile from output to input. Because of the jagged
discontinuity at RA=12 h in the HPX all-sky projection, we need to
disable wrapping of WCS values qwithint he WcsMap class temporarily. */
of the whole tile from output to input. */
lbnd_in[ 0 ] = tile_lbnd[ 0 ] - 0.5;
lbnd_in[ 1 ] = tile_lbnd[ 1 ] - 0.5;
lbnd_in[ 2 ] = lbnd[ 2 ] - 0.5;
ubnd_in[ 0 ] = tile_ubnd[ 0 ] - 0.5;
ubnd_in[ 1 ] = tile_ubnd[ 1 ] - 0.5;
ubnd_in[ 2 ] = ubnd[ 2 ] - 0.5;
old_wrap = astGetWcsWrap;
astSetWcsWrap( 0 );

astMapBox( p2pmap, lbnd_in, ubnd_in, 0, 1, lbnd_out + 0,
ubnd_out + 0, NULL, NULL );
astMapBox( p2pmap, lbnd_in, ubnd_in, 0, 2, lbnd_out + 1,
ubnd_out + 1, NULL, NULL );
if( ndim == 3 ) astMapBox( p2pmap, lbnd_in, ubnd_in, 0, 3,
lbnd_out + 2, ubnd_out + 2, NULL,
NULL );
astSetWcsWrap( old_wrap );

/* Tiles that straddle RA=12h will produce very big bounding boxes
because half of the tile will be wrapped round to the other side of the
all-sky grid. Test for this, and if found, re-calculate the bounding
boxes without wrapping WCS values within the WcsMap class. */
if( ubnd_out[0] - lbnd_out[0] > 10000 ||
ubnd_out[1] - lbnd_out[1] > 10000 ) {
old_wrap = astGetWcsWrap;
astSetWcsWrap( 0 );
astMapBox( p2pmap, lbnd_in, ubnd_in, 0, 1, lbnd_out + 0,
ubnd_out + 0, NULL, NULL );
astMapBox( p2pmap, lbnd_in, ubnd_in, 0, 2, lbnd_out + 1,
ubnd_out + 1, NULL, NULL );
if( ndim == 3 ) astMapBox( p2pmap, lbnd_in, ubnd_in, 0, 3,
lbnd_out + 2, ubnd_out + 2, NULL,
NULL );
astSetWcsWrap( old_wrap );
}

lbnd_tile[ 0 ] = floor( lbnd_out[ 0 ] ) + 1;
lbnd_tile[ 1 ] = floor( lbnd_out[ 1 ] ) + 1;
Expand Down

0 comments on commit 56bf8fa

Please sign in to comment.