Skip to content

Commit

Permalink
No Nw readding
Browse files Browse the repository at this point in the history
  • Loading branch information
@damonge
damonge committed Mar 7, 2024
1 parent c99fc55 commit 0ba15bb
Show file tree
Hide file tree
Showing 6 changed files with 36 additions and 64 deletions.
6 changes: 2 additions & 4 deletions pymaster/namaster.i
View file
Original file line number Diff line number Diff line change
Expand Up @@ -431,17 +431,15 @@ nmt_workspace *comp_coupling_matrix(int spin1,int spin2,
int nell4,double *f_ell,
nmt_binning_scheme *bin,
int is_teb,int l_toeplitz,
int l_exact,int dl_band,
double Nw)
int l_exact,int dl_band)
{
asserting(nlb1==lmax+1);
asserting(nlb2==lmax+1);
asserting(nell4==lmax_mask+1);
return nmt_compute_coupling_matrix(spin1,spin2,lmax,lmax_mask,
pure_e_1,pure_b_1,pure_e_2,pure_b_2,
f_ell,beam1,beam2,
bin,is_teb,l_toeplitz,l_exact,dl_band,
Nw);
bin,is_teb,l_toeplitz,l_exact,dl_band);
}

nmt_workspace_flat *comp_coupling_matrix_flat(nmt_field_flat *fl1,nmt_field_flat *fl2,
Expand Down
34 changes: 8 additions & 26 deletions pymaster/namaster_wrap.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3407,17 +3407,15 @@ nmt_workspace *comp_coupling_matrix(int spin1,int spin2,
int nell4,double *f_ell,
nmt_binning_scheme *bin,
int is_teb,int l_toeplitz,
int l_exact,int dl_band,
double Nw)
int l_exact,int dl_band)
{
asserting(nlb1==lmax+1);
asserting(nlb2==lmax+1);
asserting(nell4==lmax_mask+1);
return nmt_compute_coupling_matrix(spin1,spin2,lmax,lmax_mask,
pure_e_1,pure_b_1,pure_e_2,pure_b_2,
f_ell,beam1,beam2,
bin,is_teb,l_toeplitz,l_exact,dl_band,
Nw);
bin,is_teb,l_toeplitz,l_exact,dl_band);
}

nmt_workspace_flat *comp_coupling_matrix_flat(nmt_field_flat *fl1,nmt_field_flat *fl2,
Expand Down Expand Up @@ -11332,7 +11330,6 @@ SWIGINTERN PyObject *_wrap_compute_coupling_matrix(PyObject *SWIGUNUSEDPARM(self
int arg14 ;
int arg15 ;
int arg16 ;
double arg17 ;
int val1 ;
int ecode1 = 0 ;
int val2 ;
Expand Down Expand Up @@ -11365,12 +11362,10 @@ SWIGINTERN PyObject *_wrap_compute_coupling_matrix(PyObject *SWIGUNUSEDPARM(self
int ecode15 = 0 ;
int val16 ;
int ecode16 = 0 ;
double val17 ;
int ecode17 = 0 ;
PyObject *swig_obj[17] ;
PyObject *swig_obj[16] ;
nmt_workspace *result = 0 ;

if (!SWIG_Python_UnpackTuple(args, "compute_coupling_matrix", 17, 17, swig_obj)) SWIG_fail;
if (!SWIG_Python_UnpackTuple(args, "compute_coupling_matrix", 16, 16, swig_obj)) SWIG_fail;
ecode1 = SWIG_AsVal_int(swig_obj[0], &val1);
if (!SWIG_IsOK(ecode1)) {
SWIG_exception_fail(SWIG_ArgError(ecode1), "in method '" "compute_coupling_matrix" "', argument " "1"" of type '" "int""'");
Expand Down Expand Up @@ -11451,12 +11446,7 @@ SWIGINTERN PyObject *_wrap_compute_coupling_matrix(PyObject *SWIGUNUSEDPARM(self
SWIG_exception_fail(SWIG_ArgError(ecode16), "in method '" "compute_coupling_matrix" "', argument " "16"" of type '" "int""'");
}
arg16 = (int)(val16);
ecode17 = SWIG_AsVal_double(swig_obj[16], &val17);
if (!SWIG_IsOK(ecode17)) {
SWIG_exception_fail(SWIG_ArgError(ecode17), "in method '" "compute_coupling_matrix" "', argument " "17"" of type '" "double""'");
}
arg17 = (double)(val17);
result = (nmt_workspace *)nmt_compute_coupling_matrix(arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8,arg9,arg10,arg11,arg12,arg13,arg14,arg15,arg16,arg17);
result = (nmt_workspace *)nmt_compute_coupling_matrix(arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8,arg9,arg10,arg11,arg12,arg13,arg14,arg15,arg16);
resultobj = SWIG_NewPointerObj(SWIG_as_voidptr(result), SWIGTYPE_p_nmt_workspace, 0 | 0 );
return resultobj;
fail:
Expand Down Expand Up @@ -15829,7 +15819,6 @@ SWIGINTERN PyObject *_wrap_comp_coupling_matrix(PyObject *SWIGUNUSEDPARM(self),
int arg17 ;
int arg18 ;
int arg19 ;
double arg20 ;
int val1 ;
int ecode1 = 0 ;
int val2 ;
Expand Down Expand Up @@ -15862,12 +15851,10 @@ SWIGINTERN PyObject *_wrap_comp_coupling_matrix(PyObject *SWIGUNUSEDPARM(self),
int ecode18 = 0 ;
int val19 ;
int ecode19 = 0 ;
double val20 ;
int ecode20 = 0 ;
PyObject *swig_obj[17] ;
PyObject *swig_obj[16] ;
nmt_workspace *result = 0 ;

if (!SWIG_Python_UnpackTuple(args, "comp_coupling_matrix", 17, 17, swig_obj)) SWIG_fail;
if (!SWIG_Python_UnpackTuple(args, "comp_coupling_matrix", 16, 16, swig_obj)) SWIG_fail;
ecode1 = SWIG_AsVal_int(swig_obj[0], &val1);
if (!SWIG_IsOK(ecode1)) {
SWIG_exception_fail(SWIG_ArgError(ecode1), "in method '" "comp_coupling_matrix" "', argument " "1"" of type '" "int""'");
Expand Down Expand Up @@ -15969,14 +15956,9 @@ SWIGINTERN PyObject *_wrap_comp_coupling_matrix(PyObject *SWIGUNUSEDPARM(self),
SWIG_exception_fail(SWIG_ArgError(ecode19), "in method '" "comp_coupling_matrix" "', argument " "19"" of type '" "int""'");
}
arg19 = (int)(val19);
ecode20 = SWIG_AsVal_double(swig_obj[16], &val20);
if (!SWIG_IsOK(ecode20)) {
SWIG_exception_fail(SWIG_ArgError(ecode20), "in method '" "comp_coupling_matrix" "', argument " "20"" of type '" "double""'");
}
arg20 = (double)(val20);
{
try {
result = (nmt_workspace *)comp_coupling_matrix(arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8,arg9,arg10,arg11,arg12,arg13,arg14,arg15,arg16,arg17,arg18,arg19,arg20);
result = (nmt_workspace *)comp_coupling_matrix(arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8,arg9,arg10,arg11,arg12,arg13,arg14,arg15,arg16,arg17,arg18,arg19);
}
finally {
SWIG_exception(SWIG_RuntimeError,nmt_error_message);
Expand Down
8 changes: 4 additions & 4 deletions pymaster/nmtlib.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -341,8 +341,8 @@ def compute_master_coefficients(lmax, lmax_mask, npcl, pcl_masks, s1, s2, pure_e
def master_calculator_free(c):
return _nmtlib.master_calculator_free(c)

def compute_coupling_matrix(spin1, spin2, lmax, lmax_mask, pure_e1, pure_b1, pure_e2, pure_b2, pcl_masks, beam1, beam2, bin, is_teb, l_toeplitz, l_exact, dl_band, Nw):
return _nmtlib.compute_coupling_matrix(spin1, spin2, lmax, lmax_mask, pure_e1, pure_b1, pure_e2, pure_b2, pcl_masks, beam1, beam2, bin, is_teb, l_toeplitz, l_exact, dl_band, Nw)
def compute_coupling_matrix(spin1, spin2, lmax, lmax_mask, pure_e1, pure_b1, pure_e2, pure_b2, pcl_masks, beam1, beam2, bin, is_teb, l_toeplitz, l_exact, dl_band):
return _nmtlib.compute_coupling_matrix(spin1, spin2, lmax, lmax_mask, pure_e1, pure_b1, pure_e2, pure_b2, pcl_masks, beam1, beam2, bin, is_teb, l_toeplitz, l_exact, dl_band)

def update_coupling_matrix(w, n_rows, new_matrix):
return _nmtlib.update_coupling_matrix(w, n_rows, new_matrix)
Expand Down Expand Up @@ -524,8 +524,8 @@ def apomask_flat(nx, ny, lx, ly, npix_1, dout, aposize, apotype):
def synfast_new_flat(nx, ny, lx, ly, nfields, seed, ncl1, ncl2, dout):
return _nmtlib.synfast_new_flat(nx, ny, lx, ly, nfields, seed, ncl1, ncl2, dout)

def comp_coupling_matrix(spin1, spin2, lmax, lmax_mask, pure_e_1, pure_b_1, pure_e_2, pure_b_2, nlb1, nlb2, nell4, bin, is_teb, l_toeplitz, l_exact, dl_band, Nw):
return _nmtlib.comp_coupling_matrix(spin1, spin2, lmax, lmax_mask, pure_e_1, pure_b_1, pure_e_2, pure_b_2, nlb1, nlb2, nell4, bin, is_teb, l_toeplitz, l_exact, dl_band, Nw)
def comp_coupling_matrix(spin1, spin2, lmax, lmax_mask, pure_e_1, pure_b_1, pure_e_2, pure_b_2, nlb1, nlb2, nell4, bin, is_teb, l_toeplitz, l_exact, dl_band):
return _nmtlib.comp_coupling_matrix(spin1, spin2, lmax, lmax_mask, pure_e_1, pure_b_1, pure_e_2, pure_b_2, nlb1, nlb2, nell4, bin, is_teb, l_toeplitz, l_exact, dl_band)

def comp_coupling_matrix_flat(fl1, fl2, bin, lmn_x, lmx_x, lmn_y, lmx_y, is_teb):
return _nmtlib.comp_coupling_matrix_flat(fl1, fl2, bin, lmn_x, lmx_x, lmn_y, lmx_y, is_teb)
Expand Down
2 changes: 1 addition & 1 deletion pymaster/workspaces.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -137,7 +137,7 @@ def compute_coupling_matrix(self, fl1, fl2, bins, is_teb=False,
int(fl1.ainfo.lmax), int(fl1.ainfo_mask.lmax),
int(fl1.pure_e), int(fl1.pure_b), int(fl2.pure_e), int(fl2.pure_b),
fl1.beam, fl2.beam, pcl_mask.flatten()-Nw,
bins.bin, int(is_teb), l_toeplitz, l_exact, dl_band, Nw)
bins.bin, int(is_teb), l_toeplitz, l_exact, dl_band)
self.has_unbinned = True

def write_to(self, fname):
Expand Down
3 changes: 1 addition & 2 deletions src/namaster.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -745,8 +745,7 @@ nmt_workspace *nmt_compute_coupling_matrix(int spin1,int spin2,
flouble *pcl_masks,
flouble *beam1,flouble *beam2,
nmt_binning_scheme *bin,int is_teb,
int l_toeplitz,int l_exact,int dl_band,
double Nw);
int l_toeplitz,int l_exact,int dl_band);

/**
* @brief Updates the mode coupling matrix with a new one.Saves nmt_workspace structure to file
Expand Down
47 changes: 20 additions & 27 deletions src/nmt_master.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -806,8 +806,7 @@ nmt_workspace *nmt_compute_coupling_matrix(int spin1,int spin2,
flouble *pcl_masks,
flouble *beam1,flouble *beam2,
nmt_binning_scheme *bin,int is_teb,
int l_toeplitz,int l_exact,int dl_band,
double Nw)
int l_toeplitz,int l_exact,int dl_band)
{
int l2,lmax_large,lmax_fields;
nmt_workspace *w;
Expand Down Expand Up @@ -846,7 +845,7 @@ nmt_workspace *nmt_compute_coupling_matrix(int spin1,int spin2,

// Apply coupling coefficients
#pragma omp parallel default(none) \
shared(w,c,pure_e1,pure_e2,pure_b1,pure_b2,spin1,spin2,Nw)
shared(w,c,pure_e1,pure_e2,pure_b1,pure_b2,spin1,spin2)
{
int ll2,ll3;
int pe1=pure_e1,pe2=pure_e2,pb1=pure_b1,pb2=pure_b2;
Expand All @@ -858,40 +857,34 @@ nmt_workspace *nmt_compute_coupling_matrix(int spin1,int spin2,
for(ll2=0;ll2<=w->lmax;ll2++) {
for(ll3=0;ll3<=w->lmax;ll3++) {
double fac=(2*ll3+1.)*sign_overall;
double Nw_term = 0;
if((spin1 == spin2) && (ll3 >= spin1) && (ll2 >= spin1)) {
Nw_term = Nw*fac/(4*M_PI);
if(spin1 > 0)
Nw_term *= 0.5;
}
if(w->ncls==1)
w->coupling_matrix_unbinned[1*ll2+0][1*ll3+0]=fac*c->xi_00[0][ll2][ll3]+Nw_term; //TT,TT
w->coupling_matrix_unbinned[1*ll2+0][1*ll3+0]=fac*c->xi_00[0][ll2][ll3]; //TT,TT
if(w->ncls==2) {
w->coupling_matrix_unbinned[2*ll2+0][2*ll3+0]=fac*c->xi_0s[0][pe1+pe2][ll2][ll3]; //TE,TE
w->coupling_matrix_unbinned[2*ll2+1][2*ll3+1]=fac*c->xi_0s[0][pb1+pb2][ll2][ll3]; //TB,TB
}
if(w->ncls==4) {
w->coupling_matrix_unbinned[4*ll2+0][4*ll3+3]=fac*c->xi_mm[0][pe1+pe2][ll2][ll3]+Nw_term; //EE,BB
w->coupling_matrix_unbinned[4*ll2+1][4*ll3+2]=-fac*c->xi_mm[0][pe1+pb2][ll2][ll3]-Nw_term; //EB,BE
w->coupling_matrix_unbinned[4*ll2+2][4*ll3+1]=-fac*c->xi_mm[0][pb1+pe2][ll2][ll3]-Nw_term; //BE,EB
w->coupling_matrix_unbinned[4*ll2+3][4*ll3+0]=fac*c->xi_mm[0][pb1+pb2][ll2][ll3]+Nw_term; //BB,EE
w->coupling_matrix_unbinned[4*ll2+0][4*ll3+0]=fac*c->xi_pp[0][pe1+pe2][ll2][ll3]+Nw_term; //EE,EE
w->coupling_matrix_unbinned[4*ll2+1][4*ll3+1]=fac*c->xi_pp[0][pe1+pb2][ll2][ll3]+Nw_term; //EB,EB
w->coupling_matrix_unbinned[4*ll2+2][4*ll3+2]=fac*c->xi_pp[0][pb1+pe2][ll2][ll3]+Nw_term; //BE,BE
w->coupling_matrix_unbinned[4*ll2+3][4*ll3+3]=fac*c->xi_pp[0][pb1+pb2][ll2][ll3]+Nw_term; //BB,BB
w->coupling_matrix_unbinned[4*ll2+0][4*ll3+3]=fac*c->xi_mm[0][pe1+pe2][ll2][ll3]; //EE,BB
w->coupling_matrix_unbinned[4*ll2+1][4*ll3+2]=-fac*c->xi_mm[0][pe1+pb2][ll2][ll3]; //EB,BE
w->coupling_matrix_unbinned[4*ll2+2][4*ll3+1]=-fac*c->xi_mm[0][pb1+pe2][ll2][ll3]; //BE,EB
w->coupling_matrix_unbinned[4*ll2+3][4*ll3+0]=fac*c->xi_mm[0][pb1+pb2][ll2][ll3]; //BB,EE
w->coupling_matrix_unbinned[4*ll2+0][4*ll3+0]=fac*c->xi_pp[0][pe1+pe2][ll2][ll3]; //EE,EE
w->coupling_matrix_unbinned[4*ll2+1][4*ll3+1]=fac*c->xi_pp[0][pe1+pb2][ll2][ll3]; //EB,EB
w->coupling_matrix_unbinned[4*ll2+2][4*ll3+2]=fac*c->xi_pp[0][pb1+pe2][ll2][ll3]; //BE,BE
w->coupling_matrix_unbinned[4*ll2+3][4*ll3+3]=fac*c->xi_pp[0][pb1+pb2][ll2][ll3]; //BB,BB
}
if(w->ncls==7) {
w->coupling_matrix_unbinned[7*ll2+0][7*ll3+0]=fac*c->xi_00[0][ll2][ll3]+Nw_term; //TT,TT
w->coupling_matrix_unbinned[7*ll2+0][7*ll3+0]=fac*c->xi_00[0][ll2][ll3]; //TT,TT
w->coupling_matrix_unbinned[7*ll2+1][7*ll3+1]=fac*c->xi_0s[0][pe2][ll2][ll3]; //TE,TE
w->coupling_matrix_unbinned[7*ll2+2][7*ll3+2]=fac*c->xi_0s[0][pb2][ll2][ll3]; //TB,TB
w->coupling_matrix_unbinned[7*ll2+3][7*ll3+6]=fac*c->xi_mm[0][pe2+pe2][ll2][ll3]+Nw_term; //EE,BB
w->coupling_matrix_unbinned[7*ll2+4][7*ll3+5]=-fac*c->xi_mm[0][pe2+pb2][ll2][ll3]+Nw_term; //EB,BE
w->coupling_matrix_unbinned[7*ll2+5][7*ll3+4]=-fac*c->xi_mm[0][pb2+pe2][ll2][ll3]+Nw_term; //BE,EB
w->coupling_matrix_unbinned[7*ll2+6][7*ll3+3]=fac*c->xi_mm[0][pb2+pb2][ll2][ll3]+Nw_term; //BB,EE
w->coupling_matrix_unbinned[7*ll2+3][7*ll3+3]=fac*c->xi_pp[0][pe2+pe2][ll2][ll3]+Nw_term; //EE,EE
w->coupling_matrix_unbinned[7*ll2+4][7*ll3+4]=fac*c->xi_pp[0][pe2+pb2][ll2][ll3]+Nw_term; //EB,EB
w->coupling_matrix_unbinned[7*ll2+5][7*ll3+5]=fac*c->xi_pp[0][pb2+pe2][ll2][ll3]+Nw_term; //BE,BE
w->coupling_matrix_unbinned[7*ll2+6][7*ll3+6]=fac*c->xi_pp[0][pb2+pb2][ll2][ll3]+Nw_term; //BB,BB
w->coupling_matrix_unbinned[7*ll2+3][7*ll3+6]=fac*c->xi_mm[0][pe2+pe2][ll2][ll3]; //EE,BB
w->coupling_matrix_unbinned[7*ll2+4][7*ll3+5]=-fac*c->xi_mm[0][pe2+pb2][ll2][ll3]; //EB,BE
w->coupling_matrix_unbinned[7*ll2+5][7*ll3+4]=-fac*c->xi_mm[0][pb2+pe2][ll2][ll3]; //BE,EB
w->coupling_matrix_unbinned[7*ll2+6][7*ll3+3]=fac*c->xi_mm[0][pb2+pb2][ll2][ll3]; //BB,EE
w->coupling_matrix_unbinned[7*ll2+3][7*ll3+3]=fac*c->xi_pp[0][pe2+pe2][ll2][ll3]; //EE,EE
w->coupling_matrix_unbinned[7*ll2+4][7*ll3+4]=fac*c->xi_pp[0][pe2+pb2][ll2][ll3]; //EB,EB
w->coupling_matrix_unbinned[7*ll2+5][7*ll3+5]=fac*c->xi_pp[0][pb2+pe2][ll2][ll3]; //BE,BE
w->coupling_matrix_unbinned[7*ll2+6][7*ll3+6]=fac*c->xi_pp[0][pb2+pb2][ll2][ll3]; //BB,BB
}
}
} //end omp for
Expand Down

0 comments on commit 0ba15bb

Please sign in to comment.