Skip to content

Commit

Permalink
smurf: Nested indexing of JSA tiles within facets
Browse files Browse the repository at this point in the history
The tiles within each HEALPix base resolution pixel (facet)
are now numbered using the nested scheme described in the
HEALPix paper.  The offsets north-east and north-west from
the southern corner of the facet are represented by the
even (including 0th) and odd bits respectively.  This,
together with the previous commit, completes the action
item from the JSA meeting on 2013-10-28 of applying the
standard HEALPix nested numbering scheme to the JSA tiles.
  • Loading branch information
grahambell committed Oct 31, 2013
1 parent 7ed1e14 commit 5a515ff
Show file tree
Hide file tree
Showing 4 changed files with 61 additions and 18 deletions.
7 changes: 5 additions & 2 deletions applications/smurf/libsmf/smf_jsatileheader.c
Expand Up @@ -66,8 +66,11 @@
* an RA range of 9h to 15h - the top right corner of the projection
* covers the same area on the sky but has no corresponding tile). Facet
* zero contains tiles zero to N*N-1, facet one contains tiles N*N to
* 2*N*N-1, etc. Within a facet tiles are ordered raster fashion from
* bottom left to top right.
* 2*N*N-1, etc. Within a facet tiles are indexed using the "nested"
* scheme described in the HEALPix paper. This starts with pixel
* zero in the southern corner of the facet. The even bits number
* the position in the north-east direction and the odd bits number
* the position in the north-west direction.
* Authors:
* DSB: David S Berry (JAC, UCLan)
Expand Down
34 changes: 26 additions & 8 deletions applications/smurf/libsmf/smf_jsatilei2xy.c
Expand Up @@ -95,7 +95,12 @@

#include "libsmf/jsatiles.h" /* Move this to smf_typ.h and smf.h when done */


/* "Binary Magic Numbers" for selecting alternate bits. */
#define SELECT_MAGIC_0 0x55555555;
#define SELECT_MAGIC_1 0x33333333;
#define SELECT_MAGIC_2 0x0F0F0F0F;
#define SELECT_MAGIC_3 0x00FF00FF;
#define SELECT_MAGIC_4 0x0000FFFF;


void smf_jsatilei2xy( int itile, smfJSATiling *skytiling, int *xt, int *yt,
Expand Down Expand Up @@ -130,11 +135,23 @@ void smf_jsatilei2xy( int itile, smfJSATiling *skytiling, int *xt, int *yt,
/* Get the offset in tiles into this facet of the required tile. */
tj = itile - fi*nsq;

/* Within a facet, tiles are stored raster-fashion from bottom left to
top right. Get the row and column indices (zero-based) of the tile within
the facet. */
yj = tj/skytiling->ntpf;
xj = tj - yj*skytiling->ntpf;
/* Within a facet, tiles are stored in nested order from the southern
corner. Get the row and column indices (zero-based) of the tile within
the facet. xj and yj come from the even and odd bits respectively.
The bit selection method is the reverse of the "Binary Magic Numbers"
method mentioned in smf_jsatilexy2i. */
xj = tj & SELECT_MAGIC_0;
yj = (tj >> 1) & SELECT_MAGIC_0;

xj = (xj | (xj >> 1)) & SELECT_MAGIC_1;
xj = (xj | (xj >> 2)) & SELECT_MAGIC_2;
xj = (xj | (xj >> 4)) & SELECT_MAGIC_3;
xj = (xj | (xj >> 8)) & SELECT_MAGIC_4;

yj = (yj | (yj >> 1)) & SELECT_MAGIC_1;
yj = (yj | (yj >> 2)) & SELECT_MAGIC_2;
yj = (yj | (yj >> 4)) & SELECT_MAGIC_3;
yj = (yj | (yj >> 8)) & SELECT_MAGIC_4;

/* Get the offsets, in facets, along X and Y, from the bottom left tile of
the first facet to the bottom left tile of the requested facet. Note that
Expand All @@ -147,8 +164,9 @@ void smf_jsatilei2xy( int itile, smfJSATiling *skytiling, int *xt, int *yt,
fy = fc + ((11 - fi) / 8);

/* Add the offsets to the facet onto the offsets within the facet to get
the total offsets. */
*xt = fx * skytiling->ntpf + xj;
the total offsets. xj is measured north-east (left) and yj is measured
north-west (up). */
*xt = (fx + 1) * skytiling->ntpf - xj - 1;
*yt = fy * skytiling->ntpf + yj;

if (fi_) {
Expand Down
31 changes: 25 additions & 6 deletions applications/smurf/libsmf/smf_jsatilexy2i.c
Expand Up @@ -81,7 +81,11 @@

#include "libsmf/jsatiles.h" /* Move this to smf_typ.h and smf.h when done */


/* "Binary Magic Numbers" for interleaving bits. */
#define INTERLEAVE_MAGIC_1 0x00FF00FF
#define INTERLEAVE_MAGIC_2 0x0F0F0F0F
#define INTERLEAVE_MAGIC_3 0x33333333
#define INTERLEAVE_MAGIC_4 0x55555555


int smf_jsatilexy2i( int xt, int yt, smfJSATiling *skytiling, int *status ){
Expand Down Expand Up @@ -142,14 +146,29 @@ int smf_jsatilexy2i( int xt, int yt, smfJSATiling *skytiling, int *status ){
/* Get the scalar zero-based index of the first tile within the facet. */
itile = fi*skytiling->ntpf*skytiling->ntpf;

/* Get the (x,y) offsets of the tile into the facet. */
tx = xt - fx*skytiling->ntpf;
/* Get the (x,y) offsets of the tile into the facet. tx is measured
north-east (left) and ty is measured north-west (up). */
tx = (fx + 1)*skytiling->ntpf - xt - 1;
ty = yt - fy*skytiling->ntpf;

/* Get the scalar index of the tile within the facet. Add this index onto
/* Get the scalar index of the tile within the facet. Interleave the bits
of tx and ty (tx gives the even bits) and add this index onto
the index of the first tile in the facet, to get the index of the
required tile. */
itile += ty*skytiling->ntpf + tx;
required tile. Interleaving is performed using the "Binary
Magic Numbers" method with code based on the public domain
example (collection (C) 1997-2005 Sean Eron Anderson) available at:
http://graphics.stanford.edu/~seander/bithacks.html#InterleaveBMN */
tx = (tx | (tx << 8)) & INTERLEAVE_MAGIC_1;
tx = (tx | (tx << 4)) & INTERLEAVE_MAGIC_2;
tx = (tx | (tx << 2)) & INTERLEAVE_MAGIC_3;
tx = (tx | (tx << 1)) & INTERLEAVE_MAGIC_4;

ty = (ty | (ty << 8)) & INTERLEAVE_MAGIC_1;
ty = (ty | (ty << 4)) & INTERLEAVE_MAGIC_2;
ty = (ty | (ty << 2)) & INTERLEAVE_MAGIC_3;
ty = (ty | (ty << 1)) & INTERLEAVE_MAGIC_4;

itile += tx | (ty << 1);
}
}

Expand Down
7 changes: 5 additions & 2 deletions applications/smurf/libsmurf/smurf_tilelist.c
Expand Up @@ -65,8 +65,11 @@
* paper (Gorsky et. al. 2005 ApJ 622, 759) (note that the facet six is
* split equally into two triangles, one at the bottom left and one at
* the top right of the projection plane). Within a facet, tiles are
* also counted raster-fashion from bottom left to top right. All the
* tiles in the first facet come first, followed by all the tiles in
* indexed using the "nested" scheme described in the HEALPix paper.
* This starts with pixel zero in the southern corner of the facet.
* The even bits number the position in the north-east direction and
* the odd bits number the position in the north-west direction. All
* the tiles in the first facet come first, followed by all the tiles in
* the second facet, etc.
*
* This is a fairly complex scheme. To help understanding, the
Expand Down

7 comments on commit 5a515ff

@dsberry
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

But...... I thought we had made a decision to leave it as it was????

I'm worried that this opening the door to people attempting to use all sorts of external HPX software that may or may not be compatible with the way we want to do things. The idea was (or so I thought) that all tile-based operations should go through smurf commands (so that we have control of it all), and that using our own tile numbering scheme was a way to ensure this.

So have any minutes or anything been made available from this meeting on 28th October? I'd better read them....

@timj
Copy link
Member

@timj timj commented on 5a515ff Oct 31, 2013

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It was felt that using the standard tile scheme would help us with internal validation of the tiling. Makes it easier for @grahambell to do tests. JSA is only going to be using SMURF and no-one else is relevant.

@dsberry
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What is "internal validation"? What tests is Graham doing?

@timj
Copy link
Member

@timj timj commented on 5a515ff Oct 31, 2013

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Testing of the Healpix software by people within the JSA development team.

@grahambell - can I ask what testing you plan to do? The critical issues are whether the tiles abut correctly, and whether the WCS metadata associates the correct sky position with each pixel. Other things, such as equal-area pixels, heirarchical pixels, etc, seem to be very secondary concerns to me. So I'm not sure that judging the existing code by whether it is in complete agreement with another external package on such HEALPix-specific issues is that relevant. There's nothing magic about HEALPix - all we want is some way of dividing the sky up into tiles that abut easily in pixel space without any overlap, and have a describable WCS. The fact that we've adopted an HPX like system I think should be seen as a private implementation detail.

@dsberry
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

And what was the reasoning for choosing nested indexing rather than ring indexing or any other indexing scheme?

@timj
Copy link
Member

@timj timj commented on 5a515ff Oct 31, 2013

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

That one is for @grahambell to answer. I think it was based on ease of switching SMURF over from the current scheme.

@grahambell
Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The minutes are here:

https://trac.jach.hawaii.edu/trac/SCUBA2/wiki/WikiJSAmeeting

But they don't go into detail or mention ring vs. nested (those are the only two). I thought that a preference had been expressed for nested, but it did turn out that nested was also the easier change because it too is based on the facet number (*N^2) and then the index within the facet after that.

Please sign in to comment.