From 1819d58778d331b9a054744bdb32f421cace37fa Mon Sep 17 00:00:00 2001 From: Yingjie Wang Date: Sun, 7 Dec 2025 05:36:55 -0500 Subject: [PATCH] fix: if cotpsi is larger than cotpsi_max, color is black --- src/common.h | 1 + src/init.c | 1 - src/render.c | 3 +++ 3 files changed, 4 insertions(+), 1 deletion(-) diff --git a/src/common.h b/src/common.h index 4d74d22..b1b1a7d 100644 --- a/src/common.h +++ b/src/common.h @@ -9,6 +9,7 @@ #define R0 (15*Rs) #define f(r) ((1.0-Rs/((double)r))) #define bmin (1.5*sqrt(3.0)+5e-10) +#define cotpsi_max (sqrt(R0*R0/(bmin*bmin*f(R0)) - 1)) #define THETAERROR 100000 typedef struct { diff --git a/src/init.c b/src/init.c index 83bf260..c705104 100644 --- a/src/init.c +++ b/src/init.c @@ -145,7 +145,6 @@ int init(double rmax, double rela_err_limit, int *size, double **x, double **y){ double cotpsi_min = 1/rmax; sample.x[0] = cotpsi_min; sample.y[0] = chi(cotpsi_min); - double cotpsi_max = sqrt(R0*R0/(bmin*bmin*f(R0)) - 1); double chi_cotpsimax = chi(cotpsi_max); refine_interval(cotpsi_min, cotpsi_max, sample.y[0], chi_cotpsimax, rela_err_limit, &sample, 1); diff --git a/src/render.c b/src/render.c index 1ebfa03..de0d798 100644 --- a/src/render.c +++ b/src/render.c @@ -19,6 +19,7 @@ inline double xy_to_b(double x, double y) { double integrand_r(double *xy, size_t dim, void *integrand_params_void) { double cotpsi = 1.0/hypot(xy[0], xy[1]); + if (cotpsi >= cotpsi_max) return 0; Integrand_params *integrand_params = (Integrand_params *)integrand_params_void; double theta = gsl_spline_eval(integrand_params->spline_data->spline, cotpsi, integrand_params->spline_data->acc); double phi = atan2(xy[1], xy[0]); @@ -30,6 +31,7 @@ double integrand_r(double *xy, size_t dim, void *integrand_params_void) { double integrand_g(double *xy, size_t dim, void *integrand_params_void) { double cotpsi = 1.0/hypot(xy[0], xy[1]); + if (cotpsi >= cotpsi_max) return 0; Integrand_params *integrand_params = (Integrand_params *)integrand_params_void; double theta = gsl_spline_eval(integrand_params->spline_data->spline, cotpsi, integrand_params->spline_data->acc); double phi = atan2(xy[1], xy[0]); @@ -41,6 +43,7 @@ double integrand_g(double *xy, size_t dim, void *integrand_params_void) { double integrand_b(double *xy, size_t dim, void *integrand_params_void) { double cotpsi = 1.0/hypot(xy[0], xy[1]); + if (cotpsi >= cotpsi_max) return 0; Integrand_params *integrand_params = (Integrand_params *)integrand_params_void; double theta = gsl_spline_eval(integrand_params->spline_data->spline, cotpsi, integrand_params->spline_data->acc); double phi = atan2(xy[1], xy[0]);