-
Notifications
You must be signed in to change notification settings - Fork 10
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Interpolation test #232
base: master
Are you sure you want to change the base?
Interpolation test #232
Conversation
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I have left some comments that mostly explain what the code is doing at several points. I suspect that you don't want to be using this code at all for what you are doing, though. It has been extensively optimized to make the transformations between ortho-slant and HEALPix fast, but that means that it is inflexible and likely to cause you more trouble than it is worth to get it to do something different.
You would probably have better luck replacing the calls so healpix_cnv_apply
with a new function that you would write. In the new function, you could take two different obs structures, use apply_astrometry
to convert the pixel grid of the source obs structure to RA and Dec values, then use apply_astrometry
again to convert those RA/Dec values to pixel coordinates of the destination obs structure. Then you could probably use one of IDL's interpolation functions (there are two, and they work differently) to convert the values of the source image to values for the destination image.
pix2vec_ring,nside,hpx_inds,pix_coords | ||
vec2ang,pix_coords,pix_dec,pix_ra,/astro | ||
apply_astrometry, obs, ra_arr=pix_ra, dec_arr=pix_dec, x_arr=xv_hpx, y_arr=yv_hpx, /ad2xy | ||
;set up orthoslant grid | ||
x_vals = repmat(findgen(obs.dimension), 1, obs.dimension) | ||
y_vals = Transpose(x_vals) | ||
;shape them correctly | ||
x_vals = reform(x_vals, obs.dimension^2) | ||
y_vals = reform(y_vals, obs.dimension^2) | ||
;transform to ra/dec, back to orthoslant. For different observations obs would be different | ||
;between the two | ||
apply_astrometry, obs, x_arr=x_vals, y_arr=y_vals, ra_arr=ra_vals, dec_arr=dec_vals, /xy2ad | ||
apply_astrometry, obs, ra_arr=ra_vals, dec_arr=dec_vals, x_arr=xv_hpx, y_arr=yv_hpx, /ad2xy |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I don't believe you want to change the code inside this IF
statement. Restore the old code, then add your new code below.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
For manipulations of the obs
structure you should use fhd_struct_update_obs.pro
in FHD/fhd_core/setup_metadata
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
To generate consistent arrays of x and y pixel values, please use:
x_vals = meshgrid(obs.dimension, obs.elements, 1)
y_vals = meshgrid(obs.dimension, obs.elements, 2)
|
||
x_frac=1.-(xv_hpx-Floor(xv_hpx)) | ||
y_frac=1.-(yv_hpx-Floor(yv_hpx)) | ||
|
||
min_bin=Min(Floor(xv_hpx)+dimension*Floor(yv_hpx))>0L | ||
max_bin=Max(Ceil(xv_hpx)+dimension*Ceil(yv_hpx))<(dimension*elements-1L) | ||
max_bin=Max(Ceil(xv_hpx)+dimension*Ceil(yv_hpx))<(dimension*elements-1L) | ||
; give hist some data, min and max for bins, binsize=1? |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
The IDL histogram
is a very powerful tool, but is also one of the most confusing. I strongly recommend reading http://www.idlcoyote.com/tips/histogram_tutorial.html
|
||
;sa=weighting function? | ||
;set hpx inds to be the right thing for the orthoslant | ||
hpx_inds = findgen(obs.dimension^2) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Do not use findgen()
to generate array indices, ever. For larger arrays it lacks sufficient precision and you will get strange results. Use Lindgen()
instead.
|
||
;now do the FFTs back and normalize | ||
|
||
FOR fi=Long64(0),n_freq_use-1 DO *dirty_img_arr[pol_i,fi]=reform(*dirty_hpx_arr[pol_i,fi], SQRT(n_hpx), SQRT(n_hpx)) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Don't use SQRT(n_hpx)
. Instead use obs.dimension, obs.elements
@@ -151,7 +163,9 @@ FOR i=0L,n_img_use-1L DO BEGIN | |||
*ija[i]=ija0 | |||
|
|||
ENDFOR | |||
|
|||
;sa=weighting function? |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
ija
and sa
are written in the sparse matrix notation of sprsax2.pro
, in FH/fhd_core/gridding . The ija
s are the coordinates with data in the sparse matrix, and the sa
s are the values.
h00=histogram(Floor(xv_hpx)+dimension*Floor(yv_hpx),min=min_bin,max=max_bin,/binsize,reverse_ind=ri00) | ||
h01=histogram(Floor(xv_hpx)+dimension*Ceil(yv_hpx),min=min_bin,max=max_bin,/binsize,reverse_ind=ri01) | ||
h10=histogram(Ceil(xv_hpx)+dimension*Floor(yv_hpx),min=min_bin,max=max_bin,/binsize,reverse_ind=ri10) | ||
h11=histogram(Ceil(xv_hpx)+dimension*Ceil(yv_hpx),min=min_bin,max=max_bin,/binsize,reverse_ind=ri11) | ||
|
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
The above code is calculating the four corners surrounding every pixel, and putting the results in difficult-to-understand lookup tables (the "reverse indices").
No description provided.