<jupyter_start><jupyter_text>Enlib write_map bug<jupyter_code>%load_ext autoreload
%autoreload 2
from enlib import enmap,wcs as mwcs
import numpy as np
import sys,os<jupyter_output><empty_output><jupyter_text>We make a full-sky arcminute resolution geometry. I've only been able to reproduce this bug for res=1.0.<jupyter_code>res = 1.0
shape, wcs = enmap.fullsky_geometry(res=res*np.pi/180./60., proj="car")
shape = (3,)+shape<jupyter_output><empty_output><jupyter_text>We do a pix2sky that is needed by map2alm and make sure it gives a sensible result.<jupyter_code>ny, nx = shape[-2:]
vy,vx = enmap.pix2sky(shape, wcs, [np.arange(ny),np.zeros(ny)])
hy,hx = enmap.pix2sky(shape, wcs, [np.zeros(nx),np.arange(nx)])
print(vy,vx,hy,hx)<jupyter_output>(array([-1.57079633, -1.57050544, -1.57021455, ...,  1.57021455,
        1.57050544,  1.57079633]), array([-3.14159265, -3.14159265, -3.14159265, ..., -3.14159265,
       -3.14159265, -3.14159265]), array([-1.57079633, -1.57079633, -1.57079633, ..., -1.57079633,
       -1.57079633, -1.57079633]), array([-3.14159265, -3.14188354, -3.14217443, ..., -9.4239053 ,
       -9.42419618, -9.42448707]))<jupyter_text>It makes sense. We now save a map that has this geometry and load it back.<jupyter_code>root = os.environ['WORK']+"/"
enmap.write_map(root+"temp.fits",enmap.zeros(shape,wcs,dtype=np.uint8))
lshape,lwcs = enmap.read_map_geometry(root+"temp.fits")

print(shape,wcs)
print(lshape,lwcs)
print(mwcs.equal(wcs,lwcs))<jupyter_output>((3, 10801, 21600), car:{cdelt:[-0.01667,0.01667],crval:[0,0],crpix:[1.08e+04,5401]})
((3, 10801, 21600), car:{cdelt:[-0.01667,0.01667],crval:[0,0],crpix:[1.08e+04,5401]})
True<jupyter_text>The shapes and wcs of the geometry we originally made and of the saved map seem to agree. So we proceed to do the same pix2sky operation on the loaded geometry.<jupyter_code>ny, nx = lshape[-2:]
vy2,vx2 = enmap.pix2sky(lshape, lwcs, [np.arange(ny),np.zeros(ny)])
hy2,hx2 = enmap.pix2sky(lshape, lwcs, [np.zeros(nx),np.arange(nx)])
print(vy2,vx2,hy2,hx2)<jupyter_output>(array([-1.57079633, -1.57050544, -1.57021455, ...,  1.57021455,
        1.57050544,  1.57079633]), array([-3.14159265, -3.14159265, -3.14159265, ..., -3.14159265,
       -3.14159265, -3.14159265]), array([-1.57079633, -1.57079633, -1.57079633, ..., -1.57079633,
       -1.57079633, -1.57079633]), array([-3.14159265, -3.14188354, -3.14217443, ..., -9.4239053 ,
       -9.42419618, -9.42448707]))<jupyter_text>The results are all nans.<jupyter_code>print(np.all(np.isclose(vy,vy2)))
print(np.all(np.isclose(vx,vx2)))
print(np.all(np.isclose(hy,hy2)))
print(np.all(np.isclose(hx,hx2)))<jupyter_output>True
True
True
True
{'fraction_non_alphanumeric': 0.1918168168168168, 'fraction_numerical': 0.19932432432432431, 'mean_word_length': 4.935412026726058, 'pattern_counts': {'":': 0, '<': 21, '<?xml version=': 0, '>': 21, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '6', 'dir': 'jupyter-structured-clean-dedup', 'id': '492834', 'license': 'BSD-2-Clause', 'n_tokens_mistral': 1488, 'n_tokens_neox': 1181, 'n_words': 195, 'path': 'scripts/enlib-write-map-bug.ipynb', 'repo_name': 'EEmGuzman/orphics', 'size': '2664'}