From f7394ba53f68626d16cf38cee1a66d02cead9b48 Mon Sep 17 00:00:00 2001 From: olakusiak <62485729+olakusiak@users.noreply.github.com> Date: Tue, 25 Nov 2025 13:49:57 +0000 Subject: [PATCH 1/2] 2 filters for kSZ not sure this is 100% correct --- class-sz/include/class_sz.h | 6 +++ class-sz/include/class_sz_precisions.h | 1 + class-sz/include/class_sz_tools.h | 1 + class-sz/source/class_sz.c | 53 ++++++++++++++++++-------- class-sz/source/input.c | 2 + 5 files changed, 48 insertions(+), 15 deletions(-) diff --git a/class-sz/include/class_sz.h b/class-sz/include/class_sz.h index 9a5dbe86..3a3e51e7 100755 --- a/class-sz/include/class_sz.h +++ b/class-sz/include/class_sz.h @@ -1594,6 +1594,10 @@ double szcounts_ntot; double * f_unwise_filter; int unwise_filter_size; + double * l_unwise_filter2; + double * f_unwise_filter2; + int unwise_filter2_size; + double * M_min_of_z_z; double * M_min_of_z_M_min; int M_min_of_z_size; @@ -2497,6 +2501,8 @@ double get_tau_profile_at_z_m_l(double z, double get_ksz_filter_at_l(double l, struct class_sz_structure * pclass_sz); +double get_ksz_filter_at_l2(double l, + struct class_sz_structure * pclass_sz); double get_M_min_of_z(double l, struct class_sz_structure * pclass_sz); diff --git a/class-sz/include/class_sz_precisions.h b/class-sz/include/class_sz_precisions.h index 4e39785f..35efc6ca 100644 --- a/class-sz/include/class_sz_precisions.h +++ b/class-sz/include/class_sz_precisions.h @@ -38,6 +38,7 @@ class_sz_string_parameter(cib_Snu_file_nu,"/class_sz_auxiliary_files/includes/fi class_sz_string_parameter(ksz_filter_file,"/class_sz_auxiliary_files/UNWISE_galaxy_distributions/unwise_filter_functions_l_fl.txt","ksz filter file") +class_sz_string_parameter(ksz_filter_file2,"/class_sz_auxiliary_files/UNWISE_galaxy_distributions/unwise_filter_functions_l_fl2.txt","ksz filter file2") class_sz_string_parameter(ksz_template_file,"/class_sz_auxiliary_files/cl_ksz_bat.dat","ksz template file") class_sz_string_parameter(ksz_reio_template_file,"/class_sz_auxiliary_files/FBN_kSZ_PS_patchy.d.txt","ksz template file, reio contribution") diff --git a/class-sz/include/class_sz_tools.h b/class-sz/include/class_sz_tools.h index 1d89460f..fd8ea0f1 100755 --- a/class-sz/include/class_sz_tools.h +++ b/class-sz/include/class_sz_tools.h @@ -315,6 +315,7 @@ double delta_to_delta_prime_nfw( int load_unbinned_nl_yy(struct class_sz_structure * pclass_sz); int load_unbinned_nl_tt(struct class_sz_structure * pclass_sz); int load_ksz_filter(struct class_sz_structure * pclass_sz); + int load_ksz_filter2(struct class_sz_structure * pclass_sz); int load_M_min_of_z(struct class_sz_structure * pclass_sz); double get_T10_alpha_at_z(double z_asked, struct class_sz_structure * pclass_sz); double get_knl_at_z(double z_asked, struct class_sz_structure * pclass_sz); diff --git a/class-sz/source/class_sz.c b/class-sz/source/class_sz.c index 268f6cbe..5ac1c8ab 100755 --- a/class-sz/source/class_sz.c +++ b/class-sz/source/class_sz.c @@ -1133,7 +1133,7 @@ if (pclass_sz->has_kSZ_kSZ_gal_1h || pclass_sz->has_kSZ_kSZ_gal_3h || pclass_sz->has_kSZ_kSZ_gal_hf) load_ksz_filter(pclass_sz); - +load_ksz_filter2(pclass_sz); if (pclass_sz->has_tSZ_gal_1h || pclass_sz->has_tSZ_gal_2h @@ -2194,6 +2194,7 @@ double cl_tt_lensed_fft[N]; double cl_ksz_fft[N]; double cl_noise_fft[N]; double Fl_bl; +double Fl_bl2; double cl_ttf[N]; int i; double * cl_tot; @@ -2238,15 +2239,15 @@ cl_noise_fft[i] = pwl_value_1d(pclass_sz->unbinned_nl_tt_size,pclass_sz->unbinne cl_ksz_fft[i] = cl_ksz_fft[i]/l[i]/(l[i]+1.)*2.*_PI_*1e-12/pba->T_cmb/pba->T_cmb; cl_noise_fft[i] = cl_noise_fft[i]*1e-12/pba->T_cmb/pba->T_cmb; Fl_bl = get_ksz_filter_at_l(l[i],pclass_sz); - +Fl_bl2 = get_ksz_filter_at_l2(l[i],pclass_sz); double result_ttf; if (pclass_sz->compute_ksz2ksz2==1){ -result_ttf = Fl_bl*Fl_bl*(cl_ksz_fft[i])*sqrt(2./(2.*_PI_)/(2.*_PI_)); // divide by 2pi^2 and multiply by 2 -- see formula +result_ttf = Fl_bl*Fl_bl2*(cl_ksz_fft[i])*sqrt(2./(2.*_PI_)/(2.*_PI_)); // divide by 2pi^2 and multiply by 2 -- see formula } else{ -result_ttf = Fl_bl*Fl_bl*(cl_tt_lensed_fft[i]+cl_ksz_fft[i]+cl_noise_fft[i])*sqrt(2./(2.*_PI_)/(2.*_PI_)); // divide by 2pi^2 and multiply by 2 -- see formula +result_ttf = Fl_bl*Fl_bl2*(cl_tt_lensed_fft[i]+cl_ksz_fft[i]+cl_noise_fft[i])*sqrt(2./(2.*_PI_)/(2.*_PI_)); // divide by 2pi^2 and multiply by 2 -- see formula } if (l[i]ell[i]; double flprime = get_ksz_filter_at_l(lprime,pclass_sz); double abs_l_plus_lprime = sqrt(l*l+lprime*lprime+2.*l*lprime*cos(theta)); -double fl_plus_lprime = get_ksz_filter_at_l(abs_l_plus_lprime,pclass_sz); +double fl_plus_lprime = get_ksz_filter_at_l2(abs_l_plus_lprime,pclass_sz); double clprime_tt_unlensed = pwl_value_1d(ppt->l_scalar_max-2,l_unlensed,cl_tt_unlensed,lprime); if (lprimehas_kSZ_kSZ_gal_1h // free(pclass_sz->theta_kSZ2_gal_theta_grid); free(pclass_sz->l_unwise_filter); free(pclass_sz->f_unwise_filter); + free(pclass_sz->l_unwise_filter2); + free(pclass_sz->f_unwise_filter2); } if(pclass_sz->has_kSZ_kSZ_gal_1h || pclass_sz->has_kSZ_kSZ_gal_2h @@ -7959,6 +7962,7 @@ if ((int) pvectsz[pclass_sz->index_part_id_cov_hsv] == 6) { double lnk[N],lnpk[N]; int ik; double fl; + double fl2; double taul; double l; double m = exp(logM); @@ -7968,12 +7972,13 @@ if ((int) pvectsz[pclass_sz->index_part_id_cov_hsv] == 6) { lnk[ik] = log(l_min)+ik/(N-1.)*(log(l_max)-log(l_min)); l = k[ik]; fl = get_ksz_filter_at_l(l,pclass_sz); + fl2 = get_ksz_filter_at_l2(l,pclass_sz); // pvectsz[pclass_sz->index_multipole_for_tau_profile] = l; evaluate_tau_profile((l+0.5)/chi,pvecback,pvectsz,pba,pclass_sz); taul = pvectsz[pclass_sz->index_tau_profile]; Pk[ik] = fl*taul; - Pko[ik] = fl*taul; + Pko[ik] = fl2*taul; } double r[N], xi[N]; @@ -8024,6 +8029,7 @@ if ((int) pvectsz[pclass_sz->index_part_id_cov_hsv] == 6) { double lnk[N],lnpk[N]; int ik; double fl; + double fl2; double taul; double l; double m = exp(logM); @@ -8033,12 +8039,13 @@ if ((int) pvectsz[pclass_sz->index_part_id_cov_hsv] == 6) { lnk[ik] = log(l_min)+ik/(N-1.)*(log(l_max)-log(l_min)); l = k[ik]; fl = get_ksz_filter_at_l(l,pclass_sz); + fl2 = get_ksz_filter_at_l2(l,pclass_sz); // pvectsz[pclass_sz->index_multipole_for_tau_profile] = l; evaluate_tau_profile((l+0.5)/chi,pvecback,pvectsz,pba,pclass_sz); taul = pvectsz[pclass_sz->index_tau_profile]; Pk[ik] = fl*taul; - Pko[ik] = fl*taul; + Pko[ik] = fl2*taul; } double r[N], xi[N]; @@ -8088,6 +8095,7 @@ if ((int) pvectsz[pclass_sz->index_part_id_cov_hsv] == 6) { double lnk[N],lnpk[N]; int ik; double fl; + double fl2; double taul; double l; double m = exp(logM); @@ -8097,13 +8105,13 @@ if ((int) pvectsz[pclass_sz->index_part_id_cov_hsv] == 6) { lnk[ik] = log(l_min)+ik/(N-1.)*(log(l_max)-log(l_min)); l = k[ik]; fl = get_ksz_filter_at_l(l,pclass_sz); - + fl2 = get_ksz_filter_at_l2(l,pclass_sz); // pvectsz[pclass_sz->index_multipole_for_tau_profile] = l; double kp = (l+0.5)/chi; evaluate_galaxy_profile_2h(kp,m_delta_gal,r_delta_gal,c_delta_gal,pvecback,pvectsz,pba,pclass_sz); taul = pvectsz[pclass_sz->index_galaxy_profile]; Pk[ik] = fl*taul; - Pko[ik] = fl*taul; + Pko[ik] = fl2*taul; } double r[N], xi[N]; @@ -8236,6 +8244,7 @@ double k[N], Pk[N],Pko[N]; double lnk[N],lnpk[N]; int ik; double fl; +double fl2; double taul; double l; // double m = exp(logM); @@ -8245,12 +8254,12 @@ k[ik] = exp(log(l_min)+ik/(N-1.)*(log(l_max)-log(l_min))); lnk[ik] = log(l_min)+ik/(N-1.)*(log(l_max)-log(l_min)); l = k[ik]; fl = get_ksz_filter_at_l(l,pclass_sz); - +fl2 = get_ksz_filter_at_l2(l,pclass_sz); // pvectsz[pclass_sz->index_multipole_for_tau_profile] = l;//pclass_sz->ell_kSZ2_gal_multipole_grid[index_l_2]; evaluate_tau_profile((l+0.5)/chi,pvecback,pvectsz,pba,pclass_sz); taul = pvectsz[pclass_sz->index_tau_profile];//get_tau_profile_at_z_m_l(z,m,l,pclass_sz,pba); Pk[ik] = fl*taul; -Pko[ik] = fl*taul; +Pko[ik] = fl2*taul; // if(l>3e3) // printf("k = %.5e pk = %.5e\n",l,Pk[ik]); } @@ -8335,6 +8344,7 @@ double k[N], Pk[N],Pko[N]; double lnk[N],lnpk[N]; int ik; double fl; +double fl2; double taul; double l; // double m = exp(logM); @@ -8344,7 +8354,7 @@ k[ik] = exp(log(l_min)+ik/(N-1.)*(log(l_max)-log(l_min))); lnk[ik] = log(l_min)+ik/(N-1.)*(log(l_max)-log(l_min)); l = k[ik]; fl = get_ksz_filter_at_l(l,pclass_sz); - +fl2 = get_ksz_filter_at_l2(l,pclass_sz); // pvectsz[pclass_sz->index_multipole_for_tau_profile] = l;//pclass_sz->ell_kSZ2_gal_multipole_grid[index_l_2]; // evaluate_tau_profile((l+0.5)/chi,pvecback,pvectsz,pba,pclass_sz); @@ -8354,7 +8364,7 @@ fl = get_ksz_filter_at_l(l,pclass_sz); taul = pvectsz[pclass_sz->index_galaxy_profile];//get_tau_profile_at_z_m_l(z,m,l,pclass_sz,pba); Pk[ik] = fl*taul; -Pko[ik] = fl*taul; +Pko[ik] = fl2*taul; // if(l>3e3) // printf("k = %.5e pk = %.5e\n",l,Pk[ik]); } @@ -8439,6 +8449,7 @@ double k[N], Pk[N],Pko[N]; double lnk[N],lnpk[N]; int ik; double fl; +double fl2; double taul; double l; // double m = exp(logM); @@ -8448,12 +8459,12 @@ k[ik] = exp(log(l_min)+ik/(N-1.)*(log(l_max)-log(l_min))); lnk[ik] = log(l_min)+ik/(N-1.)*(log(l_max)-log(l_min)); l = k[ik]; fl = get_ksz_filter_at_l(l,pclass_sz); - +fl2 = get_ksz_filter_at_l2(l,pclass_sz); // pvectsz[pclass_sz->index_multipole_for_tau_profile] = l;//pclass_sz->ell_kSZ2_gal_multipole_grid[index_l_2]; evaluate_tau_profile((l+0.5)/chi,pvecback,pvectsz,pba,pclass_sz); taul = pvectsz[pclass_sz->index_tau_profile];//get_tau_profile_at_z_m_l(z,m,l,pclass_sz,pba); Pk[ik] = fl*taul; -Pko[ik] = fl*taul; +Pko[ik] = fl2*taul; // if(l>3e3) // printf("k = %.5e pk = %.5e\n",l,Pk[ik]); } @@ -11075,6 +11086,18 @@ else ell); return fl; } +double get_ksz_filter_at_l2(double ell, + struct class_sz_structure * pclass_sz){ + double fl = 1.; + if (ell <= pclass_sz->l_unwise_filter2[0] || ell >= pclass_sz->l_unwise_filter2[pclass_sz->unwise_filter2_size-1]) + fl = 0.; + else + fl = pwl_value_1d(pclass_sz->unwise_filter2_size, + pclass_sz->l_unwise_filter2, + pclass_sz->f_unwise_filter2, + ell); + return fl; + } double get_M_min_of_z(double ell, struct class_sz_structure * pclass_sz){ diff --git a/class-sz/source/input.c b/class-sz/source/input.c index 6e13735f..6709eec6 100755 --- a/class-sz/source/input.c +++ b/class-sz/source/input.c @@ -6557,6 +6557,8 @@ class_read_int("no_tt_noise_in_kSZ2X_cov",pclass_sz->no_tt_noise_in_kSZ2X_cov); // class_read_string("sBBN_file",ppr->sBBN_file); class_read_string("ksz_filter_file",pclass_sz->ksz_filter_file); class_read_string("projected_field_filter_file",pclass_sz->ksz_filter_file); + class_read_string("ksz_filter_file2",pclass_sz->ksz_filter_file2); + class_read_string("projected_field_filter_file2",pclass_sz->ksz_filter_file2); class_read_string("full_path_to_dndz_gal",pclass_sz->full_path_to_dndz_gal); class_read_string("full_path_and_prefix_to_dndz_ngal",pclass_sz->full_path_and_prefix_to_dndz_ngal); class_read_string("full_path_to_redshift_dependent_M_min",pclass_sz->full_path_to_redshift_dependent_M_min); From e5222e664f5be09a8b7b173e5c52b656fcb81b06 Mon Sep 17 00:00:00 2001 From: olakusiak <62485729+olakusiak@users.noreply.github.com> Date: Tue, 25 Nov 2025 13:53:59 +0000 Subject: [PATCH 2/2] forgot this file --- class-sz/tools/class_sz_tools.c | 157 +++++++++++++++++++++++++++++--- 1 file changed, 142 insertions(+), 15 deletions(-) diff --git a/class-sz/tools/class_sz_tools.c b/class-sz/tools/class_sz_tools.c index 49945f4d..0030a88d 100755 --- a/class-sz/tools/class_sz_tools.c +++ b/class-sz/tools/class_sz_tools.c @@ -13622,6 +13622,121 @@ int load_ksz_filter(struct class_sz_structure * pclass_sz) return _SUCCESS_; } +int load_ksz_filter2(struct class_sz_structure * pclass_sz) +{ + + if (pclass_sz->sz_verbose >= 1) + printf("-> loading the filter f(l) for cl^kSZ2_gal\n"); + + + class_alloc(pclass_sz->l_unwise_filter2,sizeof(double *)*100,pclass_sz->error_message); + class_alloc(pclass_sz->f_unwise_filter2,sizeof(double *)*100,pclass_sz->error_message); + //class_alloc(pclass_sz->PP_d2lnI,sizeof(double *)*100,pclass_sz->error_message); + + //char arguments[_ARGUMENT_LENGTH_MAX_]; + char line[_LINE_LENGTH_MAX_]; + //char command_with_arguments[2*_ARGUMENT_LENGTH_MAX_]; + FILE *process; + int n_data_guess, n_data = 0; + double *lnx = NULL, *lnI = NULL, *tmp = NULL; + double this_lnx, this_lnI; + int status; + int index_x; + + + /** 1. Initialization */ + /* Prepare the data (with some initial size) */ + n_data_guess = 100; + lnx = (double *)malloc(n_data_guess*sizeof(double)); + lnI = (double *)malloc(n_data_guess*sizeof(double)); + + + /** 2. Launch the command and retrieve the output */ + /* Launch the process */ + char Filepath[_ARGUMENT_LENGTH_MAX_]; + + class_open(process,pclass_sz->ksz_filter_file2, "r",pclass_sz->error_message); + if (pclass_sz->sz_verbose >= 1) + printf("-> File Name: %s\n",pclass_sz->ksz_filter_file2); + + + /* Read output and store it */ + while (fgets(line, sizeof(line)-1, process) != NULL) { + sscanf(line, "%lf %lf", &this_lnx, &this_lnI); + //printf("lnx = %e\n",this_lnx); + + + + + /* Standard technique in C: + /*if too many data, double the size of the vectors */ + /* (it is faster and safer that reallocating every new line) */ + if((n_data+1) > n_data_guess) { + n_data_guess *= 2; + tmp = (double *)realloc(lnx, n_data_guess*sizeof(double)); + class_test(tmp == NULL, + pclass_sz->error_message, + "Error allocating memory to read the pressure profile.\n"); + lnx = tmp; + tmp = (double *)realloc(lnI, n_data_guess*sizeof(double)); + class_test(tmp == NULL, + pclass_sz->error_message, + "Error allocating memory to read the pressure profile.\n"); + lnI = tmp; + }; + /* Store */ + lnx[n_data] = this_lnx; + lnI[n_data] = this_lnI; + + n_data++; + /* Check ascending order of the k's */ + if(n_data>1) { + class_test(lnx[n_data-1] <= lnx[n_data-2], + pclass_sz->error_message, + "The ell/ells's are not strictly sorted in ascending order, " + "as it is required for the calculation of the splines.\n"); + } + } + + /* Close the process */ + // status = pclose(process); + status = fclose(process); + class_test(status != 0., + pclass_sz->error_message, + "The attempt to launch the external command was unsuccessful. " + "Try doing it by hand to check for errors."); + + /** 3. Store the read results into CLASS structures */ + pclass_sz->unwise_filter2_size = n_data; + /** Make room */ + + class_realloc(pclass_sz->l_unwise_filter2, + pclass_sz->l_unwise_filter2, + pclass_sz->unwise_filter2_size*sizeof(double), + pclass_sz->error_message); + class_realloc(pclass_sz->f_unwise_filter2, + pclass_sz->f_unwise_filter2, + pclass_sz->unwise_filter2_size*sizeof(double), + pclass_sz->error_message); + + + + /** Store them */ + for (index_x=0; index_xunwise_filter2_size; index_x++) { + pclass_sz->l_unwise_filter2[index_x] = lnx[index_x]; + pclass_sz->f_unwise_filter2[index_x] = lnI[index_x]; + }; + + /** Release the memory used locally */ + free(lnx); + free(lnI); + + if (pclass_sz->sz_verbose >= 1) + printf("-> filter f(l) for cl^kSZ2_x successfully loaded.\n"); + + + return _SUCCESS_; +} @@ -18134,6 +18249,7 @@ double k[N], Pk1[N],Pk2[N], Pkr[N]; double lnk[N],lnpk[N]; int ik; double fl; +double fl2; // double taul; double l; // double m = exp(logM); @@ -18155,7 +18271,8 @@ if (isnan(Pk1[ik])||isinf(Pk1[ik])){ // pvectsz[pclass_sz->index_multipole_for_pk] = l; // evaluate_pk_at_ell_plus_one_half_over_chi(pvecback,pvectsz,pba,ppm,pnl,pclass_sz); double pkl = get_pk_lin_at_k_and_z((l+0.5)/chi,z,pba,ppm,pnl,pclass_sz);//pvectsz[pclass_sz->index_pk_for_halo_bias]; -Pk2[ik] = fl*pkl*get_psi_b1t_at_k_and_z((l+0.5)/chi,z,pclass_sz); +fl2 = get_ksz_filter_at_l2(l,pclass_sz); +Pk2[ik] = fl2*pkl*get_psi_b1t_at_k_and_z((l+0.5)/chi,z,pclass_sz); // if(l>3e3) // printf("l = %.5e pk2 = %.5e\n",l,Pk2[ik]); } @@ -18257,6 +18374,7 @@ double t12_xi12[N],t12_Pkr[N]; double lnk[N]; int ik; double fl; +double fl2; // double taul; double l; double pkl=0.; @@ -18293,6 +18411,7 @@ l = k[ik]; // pkl = pvectsz[pclass_sz->index_pk_for_halo_bias]; pkl = get_pk_lin_at_k_and_z((l+0.5)/chi,z,pba,ppm,pnl,pclass_sz); fl = get_ksz_filter_at_l(l,pclass_sz); +fl2 = get_ksz_filter_at_l2(l,pclass_sz); // if ((l+0.5)/chi>1e-2) fl = 0.; psi_bt = get_psi_b1t_at_k_and_z((l+0.5)/chi,z,pclass_sz); @@ -18304,10 +18423,10 @@ pk_phi_m2[ik] = pow((l+0.5)/chi,-2)*fl*psi_bt; pk_phi_4[ik] = pow((l+0.5)/chi,4)*fl*psi_bt; pk_phi_2[ik] = pow((l+0.5)/chi,2)*fl*psi_bt; -pk_tilde_phi_0[ik] = fl*pkl*psi_bt; -pk_tilde_phi_m2[ik] = pow((l+0.5)/chi,-2)*fl*pkl*psi_bt; -pk_tilde_phi_2[ik] = pow((l+0.5)/chi,2)*fl*pkl*psi_bt; -pk_tilde_phi_b20[ik] = fl*pkl*psi_b2t; +pk_tilde_phi_0[ik] = fl2*pkl*psi_bt; +pk_tilde_phi_m2[ik] = pow((l+0.5)/chi,-2)*fl2*pkl*psi_bt; +pk_tilde_phi_2[ik] = pow((l+0.5)/chi,2)*fl2*pkl*psi_bt; +pk_tilde_phi_b20[ik] = fl2*pkl*psi_b2t; @@ -18526,6 +18645,7 @@ double k[N], Pk1[N],Pk2[N], Pkr[N]; double lnk[N],lnpk[N]; int ik; double fl; +double fl2; // double taul; double l; // double m = exp(logM); @@ -18547,7 +18667,8 @@ if (isnan(Pk1[ik])||isinf(Pk1[ik])){ // pvectsz[pclass_sz->index_multipole_for_pk] = l; // evaluate_pk_at_ell_plus_one_half_over_chi(pvecback,pvectsz,pba,ppm,pnl,pclass_sz); double pkl = get_pk_lin_at_k_and_z((l+0.5)/chi,z,pba,ppm,pnl,pclass_sz);//pvectsz[pclass_sz->index_pk_for_halo_bias]; -Pk2[ik] = fl*pkl*get_psi_b1t_at_k_and_z((l+0.5)/chi,z,pclass_sz); +fl2 = get_ksz_filter_at_l2(l,pclass_sz); +Pk2[ik] = fl2*pkl*get_psi_b1t_at_k_and_z((l+0.5)/chi,z,pclass_sz); // if(l>3e3) // printf("k = %.5e pk = %.5e\n",l,Pk2[ik]); } @@ -18655,6 +18776,7 @@ double k[N], Pk1[N],Pk2[N], Pkr[N]; double lnk[N],lnpk[N]; int ik; double fl; +double fl2; // double taul; double l; // double m = exp(logM); @@ -18676,7 +18798,8 @@ if (isnan(Pk1[ik])||isinf(Pk1[ik])){ // pvectsz[pclass_sz->index_multipole_for_pk] = l; // evaluate_pk_at_ell_plus_one_half_over_chi(pvecback,pvectsz,pba,ppm,pnl,pclass_sz); double pkl = get_pk_lin_at_k_and_z((l+0.5)/chi,z,pba,ppm,pnl,pclass_sz);//pvectsz[pclass_sz->index_pk_for_halo_bias]; -Pk2[ik] = fl*pkl*get_psi_b1g_at_k_and_z((l+0.5)/chi,z,pclass_sz); +fl2 = get_ksz_filter_at_l2(l,pclass_sz); +Pk2[ik] = fl2*pkl*get_psi_b1g_at_k_and_z((l+0.5)/chi,z,pclass_sz); // if(l>3e3) // printf("k = %.5e pk = %.5e\n",l,Pk2[ik]); } @@ -18771,6 +18894,7 @@ double t12_xi12[N],t12_Pkr[N]; double lnk[N]; int ik; double fl; +double fl2; // double taul; double l; double pkl=0.; @@ -18818,10 +18942,11 @@ pk_phi_m2[ik] = pow((l+0.5)/chi,-2)*fl*psi_bt; pk_phi_4[ik] = pow((l+0.5)/chi,4)*fl*psi_bt; pk_phi_2[ik] = pow((l+0.5)/chi,2)*fl*psi_bt; -pk_tilde_phi_0[ik] = fl*pkl*psi_bt; -pk_tilde_phi_m2[ik] = pow((l+0.5)/chi,-2)*fl*pkl*psi_bt; -pk_tilde_phi_2[ik] = pow((l+0.5)/chi,2)*fl*pkl*psi_bt; -pk_tilde_phi_b20[ik] = fl*pkl*psi_b2t; +fl2 = get_ksz_filter_at_l2(l,pclass_sz); +pk_tilde_phi_0[ik] = fl2*pkl*psi_bt; +pk_tilde_phi_m2[ik] = pow((l+0.5)/chi,-2)*fl2*pkl*psi_bt; +pk_tilde_phi_2[ik] = pow((l+0.5)/chi,2)*fl2*pkl*psi_bt; +pk_tilde_phi_b20[ik] = fl2*pkl*psi_b2t; @@ -19028,6 +19153,7 @@ double t12_xi12[N],t12_Pkr[N]; double lnk[N]; int ik; double fl; +double fl2; // double taul; double l; double pkl=0.; @@ -19075,10 +19201,11 @@ pk_phi_m2[ik] = pow((l+0.5)/chi,-2)*fl*psi_bt; pk_phi_4[ik] = pow((l+0.5)/chi,4)*fl*psi_bt; pk_phi_2[ik] = pow((l+0.5)/chi,2)*fl*psi_bt; -pk_tilde_phi_0[ik] = fl*pkl*psi_bt; -pk_tilde_phi_m2[ik] = pow((l+0.5)/chi,-2)*fl*pkl*psi_bt; -pk_tilde_phi_2[ik] = pow((l+0.5)/chi,2)*fl*pkl*psi_bt; -pk_tilde_phi_b20[ik] = fl*pkl*psi_b2t; +fl2 = get_ksz_filter_at_l2(l,pclass_sz); +pk_tilde_phi_0[ik] = fl2*pkl*psi_bt; +pk_tilde_phi_m2[ik] = pow((l+0.5)/chi,-2)*fl2*pkl*psi_bt; +pk_tilde_phi_2[ik] = pow((l+0.5)/chi,2)*fl2*pkl*psi_bt; +pk_tilde_phi_b20[ik] = fl2*pkl*psi_b2t;