1; ----------------------------------------------------------------------------- 2; 3; Copyright (C) 1997-2013 Krzysztof M. Gorski, Eric Hivon, Anthony J. Banday 4; 5; 6; 7; 8; 9; This file is part of HEALPix. 10; 11; HEALPix is free software; you can redistribute it and/or modify 12; it under the terms of the GNU General Public License as published by 13; the Free Software Foundation; either version 2 of the License, or 14; (at your option) any later version. 15; 16; HEALPix is distributed in the hope that it will be useful, 17; but WITHOUT ANY WARRANTY; without even the implied warranty of 18; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 19; GNU General Public License for more details. 20; 21; You should have received a copy of the GNU General Public License 22; along with HEALPix; if not, write to the Free Software 23; Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA 24; 25; For more information about HEALPix see http://healpix.sourceforge.net 26; 27; ----------------------------------------------------------------------------- 28pro orth2pix, xpos, ypos, id_pix, lon_deg, lat_deg, value 29;+ 30; convert a uv position on a Orthographic map into a pixel number and 31; (long, lat) 32; only for scalar input 33; 34;- 35 36@viewcom ; define common 37do_poldirection = (polar_type eq 2) 38indlist = (n_elements(pixel_list) eq n_elements(data_plot)) 39 40if (do_flip) then flipconv = +1 else flipconv = -1 ; longitude increase leftward by default (astro convention) 41 42if (do_poldirection) then begin 43 ; compute new position of pixelisation North Pole in the plot coordinates 44 north_pole = [0.,0.,1.] 45 if (do_conv) then north_pole = SKYCONV(north_pole, inco= coord_in, outco=coord_out) 46 if (do_rot) then north_pole = north_pole # transpose(eul_mat) 47endif 48 49id_pix = -1 50lon_deg = -1000. 51lat_deg = -1000. 52 53 54if (do_fullsky) then begin 55 c0 = 1 56 c1 = -1 57 valid = ((xpos-c0)^2 + ypos^2 le 1.) or ((xpos-c1)^2 + ypos^2 le 1.) 58endif else begin 59 valid = (xpos)^2 + ypos^2 le 1. 60endelse 61 62if (valid) then begin 63 64 nested = keyword_set(nest) 65 v1 = ypos 66 u1 = flipconv * xpos 67 if (do_fullsky) then begin 68 sign = (u1 ge 0.)*2 - 1 ; =1 for u1>0 , =-1 otherwise 69 ys = abs(u1)-c0 70 xs = sqrt(1.-ys*ys-v1*v1) 71 vector = [[sign*xs],[ys],[v1]] ; normalized vector 72 endif else begin 73 xs = sqrt(1.-u1*u1-v1*v1) 74 vector = [[xs],[u1],[v1]] ; normalized vector 75 endelse 76 ; --------- 77 ; rotation 78 ; --------- 79 if (do_rot) then vector = vector # eul_mat 80 vec2ang, vector, lat_deg, lon_deg, /astro 81 if (do_conv) then vector = SKYCONV(vector, inco = coord_out, outco = coord_in) 82 ; we go from the final Orthographic map (system coord_out) to 83 ; the original one (system coord_in) 84 ; ------------------------------------------------------------- 85 ; converts the position on the sphere into pixel number 86 ; and project the corresponding data value on the map 87 ; ------------------------------------------------------------- 88 case pix_type of 89 'R' : begin 90 vec2pix_ring, pix_param, vector, id_pix ; Healpix ring 91 pix2vec_ring, pix_param, id_pix, vectorp ; vector pointing to Healpix pixel center 92 end 93 'N' : begin 94 vec2pix_nest, pix_param, vector, id_pix ; Healpix nest 95 pix2vec_nest, pix_param, id_pix, vectorp ; vector pointing to Healpix pixel center 96 end 97 'Q' : begin 98 id_pix = uv2pix(vector, pix_param) ; QuadCube (COBE cgis software) 99 vectorp = pix2uv(id_pix, pix_param) ; vector pointing to QuadCube pixel center 100 end 101 else : print,'error on pix_type' 102 endcase 103 ; -------------------------------- 104 ; deal with polarisation direction 105 ; -------------------------------- 106 if (do_poldirection) then begin 107 phi = 0. 108 if (do_rot or do_conv) then begin 109 ; rotate pixel center to final coordinate system 110 if (do_conv) then vectorp = SKYCONV(vectorp, inco= coord_in, outco=coord_out) 111 if (do_rot) then vectorp = vectorp # transpose(eul_mat) 112 ; compute rotation of local coordinates around each vector 113 tmp_sin = north_pole[1] * vectorp[*,0] - north_pole[0] * vectorp[*,1] 114 tmp_cos = north_pole[2] - vectorp[*,2] * (north_pole[0:2] ## vectorp) 115 if (flipconv lt 0) then tmp_cos = flipconv * tmp_cos 116 phi = ATAN(tmp_sin, tmp_cos) ; angle in radians 117 endif 118 value = (data_plot[id_pix] - phi + 4*!PI) MOD (2*!PI) ; in 0,2pi 119 endif else begin 120 if (indlist) then begin 121 value = sample_sparse_array(data_plot, id_pix, in=pixel_list, default=!healpix.bad_value) 122 endif else begin 123 value = data_plot[id_pix] 124 endelse 125 endelse 126 if (n_elements(xpos) eq 1) then begin 127 id_pix = id_pix[0] 128 value = value[0] 129 endif 130endif 131 132return 133end 134