s2let_hpx_plot_mollweide

 s2let_hpx_plot_mollweide 
 Plot a real Healpix map using Mollweide projection.

 S2LET package to perform Wavelets transform on the Sphere.
 Copyright (C) 2012-2015  Boris Leistedt & Jason McEwen
 See LICENSE.txt for license details

0001 function s2let_hpx_plot_mollweide(f)
0002 
0003 % s2let_hpx_plot_mollweide
0004 % Plot a real Healpix map using Mollweide projection.
0005 %
0006 % S2LET package to perform Wavelets transform on the Sphere.
0007 % Copyright (C) 2012-2015  Boris Leistedt & Jason McEwen
0008 % See LICENSE.txt for license details
0009 
0010 sz = size(f);
0011 nside = sqrt(max(sz)/12);
0012 
0013 [thetas, phis] = s2let_hpx_sampling_ring(nside);
0014           
0015 [x, y] = ssht_mollweide(thetas, phis);
0016 
0017 gridDelaunay = delaunay(x,y);
0018 h = trisurf(gridDelaunay,x,y,f*0.0,f);
0019 
0020 set(h, 'LineStyle', 'none')
0021 axis equal
0022 axis off
0023 campos([0 0 1])
0024 camup([0 1 0])
0025 
0026 
0027 function [x, y] = ssht_mollweide(thetas, phis)
0028 % ssht_mollweide - Compute Mollweide projection
0029 %
0030 % Compute Mollweide projection of spherical coordinates.
0031 %
0032 % Usage is given by
0033 %
0034 %   [x,y] = ssht_mollweide(thetas, phis)
0035 %
0036 % where thetas and phis are spherical coordinates and x and y are the
0037 % projected Mollweide coordinates.
0038 %
0039 % Author: Jason McEwen (www.jasonmcewen.org)
0040 
0041 MAX_ITERATIONS = 1e5;
0042 TOL = 1e-10;
0043 
0044 % Convert theta to longitude.
0045 thetas = pi/2 - thetas;
0046 phis = phis - pi;
0047 
0048 t = thetas;
0049 for it = 1:MAX_ITERATIONS
0050 
0051    dt = (t + sin(t) - pi.*sin(thetas)) ./ (1 + cos(t));
0052    t = t - dt;
0053    
0054    if(max(abs(dt)) < TOL)
0055       break;
0056    end
0057    
0058 end
0059 t = t/2;
0060 x = 2 .* sqrt(2) ./ pi .* phis .* cos(t);
0061 y = sqrt(2) .* sin(t);