From 3b408b4b00794090d4c04f1d7fb95b4e7700a835 Mon Sep 17 00:00:00 2001 From: lidihei <912219465@qq.com> Date: Wed, 19 Oct 2022 14:30:16 +0800 Subject: [PATCH 01/12] add model from arguments in trm_lcurve.cc --- src/trm_lcurve.cc | 267 ++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 267 insertions(+) diff --git a/src/trm_lcurve.cc b/src/trm_lcurve.cc index 3887d77..3cf5bfd 100644 --- a/src/trm_lcurve.cc +++ b/src/trm_lcurve.cc @@ -417,6 +417,273 @@ Lcurve::Model::Model(const std::string& file) { } + + + +/* Constructor from parameters added by lijiao + */ +Lcurve::Model::Model(//Binary and stars + double q_value, double q_range, double q_dstep, bool q_vary, bool q_defined, + double iangle_value, double iangle_range, double iangle_dstep, bool iangle_vary, bool iangle_defined, + double r1_value, double r1_range, double r1_dstep, bool r1_vary, bool r1_defined, + double r2_value, double r2_range, double r2_dstep, bool r2_vary, bool r2_defined, + double cphi3_value, double cphi3_range, double cphi3_dstep, bool cphi3_vary, bool cphi3_defined, + double cphi4_value, double cphi4_range, double cphi4_dstep, bool cphi4_vary, bool cphi4_defined, + double spin1_value, double spin1_range, double spin1_dstep, bool spin1_vary, bool spin1_defined, + double spin2_value, double spin2_range, double spin2_dstep, bool spin2_vary, bool spin2_defined, + double t1_value, double t1_range, double t1_dstep, bool t1_vary, bool t1_defined, + double t2_value, double t2_range, double t2_dstep, bool t2_vary, bool t2_defined, + double ldc1_1_value, double ldc1_1_range, double ldc1_1_dstep, bool ldc1_1_vary, bool ldc1_1_defined, + double ldc1_2_value, double ldc1_2_range, double ldc1_2_dstep, bool ldc1_2_vary, bool ldc1_2_defined, + double ldc1_3_value, double ldc1_3_range, double ldc1_3_dstep, bool ldc1_3_vary, bool ldc1_3_defined, + double ldc1_4_value, double ldc1_4_range, double ldc1_4_dstep, bool ldc1_4_vary, bool ldc1_4_defined, + double ldc2_1_value, double ldc2_1_range, double ldc2_1_dstep, bool ldc2_1_vary, bool ldc2_1_defined, + double ldc2_2_value, double ldc2_2_range, double ldc2_2_dstep, bool ldc2_2_vary, bool ldc2_2_defined, + double ldc2_3_value, double ldc2_3_range, double ldc2_3_dstep, bool ldc2_3_vary, bool ldc2_3_defined, + double ldc2_4_value, double ldc2_4_range, double ldc2_4_dstep, bool ldc2_4_vary, bool ldc2_4_defined, + double velocity_scale_value, double velocity_scale_range, double velocity_scale_dstep, bool velocity_scale_vary, bool velocity_scale_defined, + double beam_factor1_value, double beam_factor1_range, double beam_factor1_dstep, bool beam_factor1_vary, bool beam_factor1_defined, + double beam_factor2_value, double beam_factor2_range, double beam_factor2_dstep, bool beam_factor2_vary, bool beam_factor2_defined, + //General + double t0_value, double t0_range, double t0_dstep, bool t0_vary, bool t0_defined, + double period_value, double period_range, double period_dstep, bool period_vary, bool period_defined, + double pdot_value, double pdot_range, double pdot_dstep, bool pdot_vary, bool pdot_defined, + double deltat_value, double deltat_range, double deltat_dstep, bool deltat_vary, bool deltat_defined, + double gravity_dark1_value, double gravity_dark1_range, double gravity_dark1_dstep, bool gravity_dark1_vary, bool gravity_dark1_defined, + double gravity_dark2_value, double gravity_dark2_range, double gravity_dark2_dstep, bool gravity_dark2_vary, bool gravity_dark2_defined, + double absorb_value, double absorb_range, double absorb_dstep, bool absorb_vary, bool absorb_defined, + double slope_value, double slope_range, double slope_dstep, bool slope_vary, bool slope_defined, + double quad_value, double quad_range, double quad_dstep, bool quad_vary, bool quad_defined, + double cube_value, double cube_range, double cube_dstep, bool cube_vary, bool cube_defined, + double third_value, double third_range, double third_dstep, bool third_vary, bool third_defined, + // Star spots + double stsp11_long_value, double stsp11_long_range, double stsp11_long_dstep, bool stsp11_long_vary, bool stsp11_long_defined, + double stsp11_lat_value, double stsp11_lat_range, double stsp11_lat_dstep, bool stsp11_lat_vary, bool stsp11_lat_defined, + double stsp11_fwhm_value, double stsp11_fwhm_range, double stsp11_fwhm_dstep, bool stsp11_fwhm_vary, bool stsp11_fwhm_defined, + double stsp11_tcen_value, double stsp11_tcen_range, double stsp11_tcen_dstep, bool stsp11_tcen_vary, bool stsp11_tcen_defined, + double stsp12_long_value, double stsp12_long_range, double stsp12_long_dstep, bool stsp12_long_vary, bool stsp12_long_defined, + double stsp12_lat_value, double stsp12_lat_range, double stsp12_lat_dstep, bool stsp12_lat_vary, bool stsp12_lat_defined, + double stsp12_fwhm_value, double stsp12_fwhm_range, double stsp12_fwhm_dstep, bool stsp12_fwhm_vary, bool stsp12_fwhm_defined, + double stsp12_tcen_value, double stsp12_tcen_range, double stsp12_tcen_dstep, bool stsp12_tcen_vary, bool stsp12_tcen_defined, + double stsp13_long_value, double stsp13_long_range, double stsp13_long_dstep, bool stsp13_long_vary, bool stsp13_long_defined, + double stsp13_lat_value, double stsp13_lat_range, double stsp13_lat_dstep, bool stsp13_lat_vary, bool stsp13_lat_defined, + double stsp13_fwhm_value, double stsp13_fwhm_range, double stsp13_fwhm_dstep, bool stsp13_fwhm_vary, bool stsp13_fwhm_defined, + double stsp13_tcen_value, double stsp13_tcen_range, double stsp13_tcen_dstep, bool stsp13_tcen_vary, bool stsp13_tcen_defined, + double stsp21_long_value, double stsp21_long_range, double stsp21_long_dstep, bool stsp21_long_vary, bool stsp21_long_defined, + double stsp21_lat_value, double stsp21_lat_range, double stsp21_lat_dstep, bool stsp21_lat_vary, bool stsp21_lat_defined, + double stsp21_fwhm_value, double stsp21_fwhm_range, double stsp21_fwhm_dstep, bool stsp21_fwhm_vary, bool stsp21_fwhm_defined, + double stsp21_tcen_value, double stsp21_tcen_range, double stsp21_tcen_dstep, bool stsp21_tcen_vary, bool stsp21_tcen_defined, + double stsp22_long_value, double stsp22_long_range, double stsp22_long_dstep, bool stsp22_long_vary, bool stsp22_long_defined, + double stsp22_lat_value, double stsp22_lat_range, double stsp22_lat_dstep, bool stsp22_lat_vary, bool stsp22_lat_defined, + double stsp22_fwhm_value, double stsp22_fwhm_range, double stsp22_fwhm_dstep, bool stsp22_fwhm_vary, bool stsp22_fwhm_defined, + double stsp22_tcen_value, double stsp22_tcen_range, double stsp22_tcen_dstep, bool stsp22_tcen_vary, bool stsp22_tcen_defined, + double uesp_long1_value, double uesp_long1_range, double uesp_long1_dstep, bool uesp_long1_vary, bool uesp_long1_defined, + double uesp_long2_value, double uesp_long2_range, double uesp_long2_dstep, bool uesp_long2_vary, bool uesp_long2_defined, + double uesp_lathw_value, double uesp_lathw_range, double uesp_lathw_dstep, bool uesp_lathw_vary, bool uesp_lathw_defined, + double uesp_taper_value, double uesp_taper_range, double uesp_taper_dstep, bool uesp_taper_vary, bool uesp_taper_defined, + double uesp_temp_value, double uesp_temp_range, double uesp_temp_dstep, bool uesp_temp_vary, bool uesp_temp_defined, + // disc + double rdisc1_value, double rdisc1_range, double rdisc1_dstep, bool rdisc1_vary, bool rdisc1_defined, //disc + double rdisc2_value, double rdisc2_range, double rdisc2_dstep, bool rdisc2_vary, bool rdisc2_defined, + double height_disc_value, double height_disc_range, double height_disc_dstep, bool height_disc_vary, bool height_disc_defined, + double beta_disc_value, double beta_disc_range, double beta_disc_dstep, bool beta_disc_vary, bool beta_disc_defined, + double temp_disc_value, double temp_disc_range, double temp_disc_dstep, bool temp_disc_vary, bool temp_disc_defined, + double texp_disc_value, double texp_disc_range, double texp_disc_dstep, bool texp_disc_vary, bool texp_disc_defined, + double lin_limb_disc_value, double lin_limb_disc_range, double lin_limb_disc_dstep, bool lin_limb_disc_vary, bool lin_limb_disc_defined, + double quad_limb_disc_value, double quad_limb_disc_range, double quad_limb_disc_dstep, bool quad_limb_disc_vary, bool quad_limb_disc_defined, + double temp_edge_value, double temp_edge_range, double temp_edge_dstep, bool temp_edge_vary, bool temp_edge_defined, + double absorb_edge_value, double absorb_edge_range, double absorb_edge_dstep, bool absorb_edge_vary, bool absorb_edge_defined, + //Bright-spot + double radius_spot_value, double radius_spot_range, double radius_spot_dstep, bool radius_spot_vary, bool radius_spot_defined, //Bright-spot + double length_spot_value, double length_spot_range, double length_spot_dstep, bool length_spot_vary, bool length_spot_defined, + double height_spot_value, double height_spot_range, double height_spot_dstep, bool height_spot_vary, bool height_spot_defined, + double expon_spot_value, double expon_spot_range, double expon_spot_dstep, bool expon_spot_vary, bool expon_spot_defined, + double epow_spot_value, double epow_spot_range, double epow_spot_dstep, bool epow_spot_vary, bool epow_spot_defined, + double angle_spot_value, double angle_spot_range, double angle_spot_dstep, bool angle_spot_vary, bool angle_spot_defined, + double yaw_spot_value, double yaw_spot_range, double yaw_spot_dstep, bool yaw_spot_vary, bool yaw_spot_defined, + double temp_spot_value, double temp_spot_range, double temp_spot_dstep, bool temp_spot_vary, bool temp_spot_defined, + double tilt_spot_value, double tilt_spot_range, double tilt_spot_dstep, bool tilt_spot_vary, bool tilt_spot_defined, + double cfrac_spot_value, double cfrac_spot_range, double cfrac_spot_dstep, bool cfrac_spot_vary, bool cfrac_spot_defined, + // Computational parameters + double delta_phase, int nlat1f, int nlat2f, int nlat1c, int nlat2c, bool npole, + int nlatfill, int nlngfill, double lfudge, double llo, double lhi, double phase1, double phase2, int nrad, double wavelength, + bool roche1, bool roche2, bool eclipse1, bool eclipse2, bool glens1, bool use_radii, + double tperiod, double gdark_bolom1, double gdark_bolom2, double mucrit1, double mucrit2, + const char* pslimb1, const char* pslimb2, bool mirror, bool add_disc, bool opaque, bool add_spot, int nspot, bool iscale + ) { + + + // Star spots: + // + // 11 = star 1, spot 1 + // 21 = star 2, spot 1 + // + + // Computational + + // Parameter value pairs + // Initialise physical parameters + q = Pparam(q_value, q_range, q_dstep, q_vary, q_defined); + // std::cout << "lijiao trm_lroche.cc: q2): " << q.value << std::endl; + iangle = Pparam(iangle_value, iangle_range, iangle_dstep, iangle_vary, iangle_defined); + r1 = Pparam(r1_value, r1_range, r1_dstep, r1_vary, r1_defined); + r2 = Pparam(r2_value, r2_range, r2_dstep, r2_vary, r2_defined); + cphi3 = Pparam(cphi3_value, cphi3_range, cphi3_dstep, cphi3_vary, cphi3_defined); + cphi4 = Pparam(cphi4_value, cphi4_range, cphi4_dstep, cphi4_vary, cphi4_defined); + spin1 = Pparam(spin1_value, spin1_range, spin1_dstep, spin1_vary, spin1_defined); + spin2 = Pparam(spin2_value, spin2_range, spin2_dstep, spin2_vary, spin2_defined); + t1 = Pparam(t1_value, t1_range, t1_dstep, t1_vary, t1_defined); + t2 = Pparam(t2_value, t2_range, t2_dstep, t2_vary, t2_defined); + ldc1_1 = Pparam(ldc1_1_value, ldc1_1_range, ldc1_1_dstep, ldc1_1_vary, ldc1_1_defined); + ldc1_2 = Pparam(ldc1_2_value, ldc1_2_range, ldc1_2_dstep, ldc1_2_vary, ldc1_2_defined); + ldc1_3 = Pparam(ldc1_3_value, ldc1_3_range, ldc1_3_dstep, ldc1_3_vary, ldc1_3_defined); + ldc1_4 = Pparam(ldc1_4_value, ldc1_4_range, ldc1_4_dstep, ldc1_4_vary, ldc1_4_defined); + ldc2_1 = Pparam(ldc2_1_value, ldc2_1_range, ldc2_1_dstep, ldc2_1_vary, ldc2_1_defined); + ldc2_2 = Pparam(ldc2_2_value, ldc2_2_range, ldc2_2_dstep, ldc2_2_vary, ldc2_2_defined); + ldc2_3 = Pparam(ldc2_3_value, ldc2_3_range, ldc2_3_dstep, ldc2_3_vary, ldc2_3_defined); + ldc2_4 = Pparam(ldc2_4_value, ldc2_4_range, ldc2_4_dstep, ldc2_4_vary, ldc2_4_defined); + velocity_scale = Pparam(velocity_scale_value, velocity_scale_range, velocity_scale_dstep, velocity_scale_vary, velocity_scale_defined); + beam_factor1 = Pparam(beam_factor1_value, beam_factor1_range, beam_factor1_dstep, beam_factor1_vary, beam_factor1_defined); + beam_factor2 = Pparam(beam_factor2_value, beam_factor2_range, beam_factor2_dstep, beam_factor2_vary, beam_factor2_defined); + t0 = Pparam(t0_value, t0_range, t0_dstep, t0_vary, t0_defined); + period = Pparam(period_value, period_range, period_dstep, period_vary, period_defined); + pdot = Pparam(pdot_value, pdot_range, pdot_dstep, pdot_vary, pdot_defined); + deltat = Pparam(deltat_value, deltat_range, deltat_dstep, deltat_vary, deltat_defined); + gravity_dark1 = Pparam(gravity_dark1_value, gravity_dark1_range, gravity_dark1_dstep, gravity_dark1_vary, gravity_dark1_defined); + gravity_dark2 = Pparam(gravity_dark2_value, gravity_dark2_range, gravity_dark2_dstep, gravity_dark2_vary, gravity_dark2_defined); + absorb = Pparam(absorb_value, absorb_range, absorb_dstep, absorb_vary, absorb_defined); + slope = Pparam(slope_value, slope_range, slope_dstep, slope_vary, slope_defined); + quad = Pparam(quad_value, quad_range, quad_dstep, quad_vary, quad_defined); + cube = Pparam(cube_value, cube_range, cube_dstep, cube_vary, cube_defined); + third = Pparam(third_value, third_range, third_dstep, third_vary, third_defined); + rdisc1 = Pparam(rdisc1_value, rdisc1_range, rdisc1_dstep, rdisc1_vary, rdisc1_defined); + rdisc2 = Pparam(rdisc2_value, rdisc2_range, rdisc2_dstep, rdisc2_vary, rdisc2_defined); + height_disc = Pparam(height_disc_value, height_disc_range, height_disc_dstep, height_disc_vary, height_disc_defined); + beta_disc = Pparam(beta_disc_value, beta_disc_range, beta_disc_dstep, beta_disc_vary, beta_disc_defined); + temp_disc = Pparam(temp_disc_value, temp_disc_range, temp_disc_dstep, temp_disc_vary, temp_disc_defined); + texp_disc = Pparam(texp_disc_value, texp_disc_range, texp_disc_dstep, texp_disc_vary, texp_disc_defined); + lin_limb_disc = Pparam(lin_limb_disc_value, lin_limb_disc_range, lin_limb_disc_dstep, lin_limb_disc_vary, lin_limb_disc_defined); + quad_limb_disc = Pparam(quad_limb_disc_value, quad_limb_disc_range, quad_limb_disc_dstep, quad_limb_disc_vary, quad_limb_disc_defined); + temp_edge = Pparam(temp_edge_value, temp_edge_range, temp_edge_dstep, temp_edge_vary, temp_edge_defined); + absorb_edge = Pparam(absorb_edge_value, absorb_edge_range, absorb_edge_dstep, absorb_edge_vary, absorb_edge_defined); + radius_spot = Pparam(radius_spot_value, radius_spot_range, radius_spot_dstep, radius_spot_vary, radius_spot_defined); + length_spot = Pparam(length_spot_value, length_spot_range, length_spot_dstep, length_spot_vary, length_spot_defined); + height_spot = Pparam(height_spot_value, height_spot_range, height_spot_dstep, height_spot_vary, height_spot_defined); + expon_spot = Pparam(expon_spot_value, expon_spot_range, expon_spot_dstep, expon_spot_vary, expon_spot_defined); + epow_spot = Pparam(epow_spot_value, epow_spot_range, epow_spot_dstep, epow_spot_vary, epow_spot_defined); + angle_spot = Pparam(angle_spot_value, angle_spot_range, angle_spot_dstep, angle_spot_vary, angle_spot_defined); + yaw_spot = Pparam(yaw_spot_value, yaw_spot_range, yaw_spot_dstep, yaw_spot_vary, yaw_spot_defined); + temp_spot = Pparam(temp_spot_value, temp_spot_range, temp_spot_dstep, temp_spot_vary, temp_spot_defined); + tilt_spot = Pparam(tilt_spot_value, tilt_spot_range, tilt_spot_dstep, tilt_spot_vary, tilt_spot_defined); + cfrac_spot = Pparam(cfrac_spot_value, cfrac_spot_range, cfrac_spot_dstep, cfrac_spot_vary, cfrac_spot_defined); + + // star-spot parameters need not have been defined, but + // if one of a group has, then all of them should be + if(stsp11_long_defined || stsp11_lat_defined || + stsp11_fwhm_defined || stsp11_tcen_defined){ + if(!(stsp11_long_defined && stsp11_lat_defined && + stsp11_fwhm_defined && stsp11_tcen_defined)) + throw Lcurve_Error("One or more of the star spot 11 parameters were not initialised"); + stsp11_long = Pparam(stsp11_long_value, stsp11_long_range, stsp11_long_dstep, stsp11_long_vary, stsp11_long_defined); + stsp11_lat = Pparam(stsp11_lat_value, stsp11_lat_range, stsp11_lat_dstep, stsp11_lat_vary, stsp11_lat_defined); + stsp11_fwhm = Pparam(stsp11_fwhm_value, stsp11_fwhm_range, stsp11_fwhm_dstep, stsp11_fwhm_vary, stsp11_fwhm_defined); + stsp11_tcen = Pparam(stsp11_tcen_value, stsp11_tcen_range, stsp11_tcen_dstep, stsp11_tcen_vary, stsp11_tcen_defined); + } + + if(stsp12_long_defined || stsp12_lat_defined || + stsp12_fwhm_defined || stsp12_tcen_defined){ + if(!(stsp12_long_defined && stsp12_lat_defined && + stsp12_fwhm_defined && stsp12_tcen_defined)) + throw Lcurve_Error("One or more of the star spot 12 parameters were not initialised"); + stsp12_long = Pparam(stsp12_long_value, stsp12_long_range, stsp12_long_dstep, stsp12_long_vary, stsp12_long_defined); + stsp12_lat = Pparam(stsp12_lat_value, stsp12_lat_range, stsp12_lat_dstep, stsp12_lat_vary, stsp12_lat_defined); + stsp12_fwhm = Pparam(stsp12_fwhm_value, stsp12_fwhm_range, stsp12_fwhm_dstep, stsp12_fwhm_vary, stsp12_fwhm_defined); + stsp12_tcen = Pparam(stsp12_tcen_value, stsp12_tcen_range, stsp12_tcen_dstep, stsp12_tcen_vary, stsp12_tcen_defined); + } + + if(stsp13_long_defined || stsp13_lat_defined || + stsp13_fwhm_defined || stsp13_tcen_defined){ + if(!(stsp13_long_defined && stsp13_lat_defined && + stsp13_fwhm_defined && stsp13_tcen_defined)) + throw Lcurve_Error("One or more of the star spot 13 parameters were not initialised"); + stsp13_long = Pparam(stsp13_long_value, stsp13_long_range, stsp13_long_dstep, stsp13_long_vary, stsp13_long_defined); + stsp13_lat = Pparam(stsp13_lat_value, stsp13_lat_range, stsp13_lat_dstep, stsp13_lat_vary, stsp13_lat_defined); + stsp13_fwhm = Pparam(stsp13_fwhm_value, stsp13_fwhm_range, stsp13_fwhm_dstep, stsp13_fwhm_vary, stsp13_fwhm_defined); + stsp13_tcen = Pparam(stsp13_tcen_value, stsp13_tcen_range, stsp13_tcen_dstep, stsp13_tcen_vary, stsp13_tcen_defined); + } + + if(stsp21_long_defined || stsp21_lat_defined || + stsp21_fwhm_defined || stsp21_tcen_defined){ + if(!(stsp21_long_defined && stsp21_lat_defined && + stsp21_fwhm_defined && stsp21_tcen_defined)) + throw Lcurve_Error("One or more of the star spot 21 parameters were not initialised"); + stsp21_long = Pparam(stsp21_long_value, stsp21_long_range, stsp21_long_dstep, stsp21_long_vary, stsp21_long_defined); + stsp21_lat = Pparam(stsp21_lat_value, stsp21_lat_range, stsp21_lat_dstep, stsp21_lat_vary, stsp21_lat_defined); + stsp21_fwhm = Pparam(stsp21_fwhm_value, stsp21_fwhm_range, stsp21_fwhm_dstep, stsp21_fwhm_vary, stsp21_fwhm_defined); + stsp21_tcen = Pparam(stsp21_tcen_value, stsp21_tcen_range, stsp21_tcen_dstep, stsp21_tcen_vary, stsp21_tcen_defined); + } + + if(stsp22_long_defined || stsp22_lat_defined || + stsp22_fwhm_defined || stsp22_tcen_defined){ + if(!(stsp22_long_defined && stsp22_lat_defined && + stsp22_fwhm_defined && stsp22_tcen_defined)) + throw Lcurve_Error("One or more of the star spot 22 parameters were not initialised"); + stsp22_long = Pparam(stsp22_long_value, stsp22_long_range, stsp22_long_dstep, stsp22_long_vary, stsp22_long_defined); + stsp22_lat = Pparam(stsp22_lat_value, stsp22_lat_range, stsp22_lat_dstep, stsp22_lat_vary, stsp22_lat_defined); + stsp22_fwhm = Pparam(stsp22_fwhm_value, stsp22_fwhm_range, stsp22_fwhm_dstep, stsp22_fwhm_vary, stsp22_fwhm_defined); + stsp22_tcen = Pparam(stsp22_tcen_value, stsp22_tcen_range, stsp22_tcen_dstep, stsp22_tcen_vary, stsp22_tcen_defined); + } + + if(uesp_long1_defined || uesp_long2_defined || + uesp_lathw_defined || uesp_taper_defined || uesp_temp_defined){ + if(!(uesp_long1_defined && uesp_long2_defined && + uesp_lathw_defined && uesp_taper_defined && uesp_temp_defined)) + throw Lcurve_Error("One or more of uniform equatorial spot parameters were not initialised"); + uesp_long1 = Pparam(uesp_long1_value, uesp_long1_range, uesp_long1_dstep, uesp_long1_vary, uesp_long1_defined); + uesp_long2 = Pparam(uesp_long2_value, uesp_long2_range, uesp_long2_dstep, uesp_long2_vary, uesp_long2_defined); + uesp_lathw = Pparam(uesp_lathw_value, uesp_lathw_range, uesp_lathw_dstep, uesp_lathw_vary, uesp_lathw_defined); + uesp_taper = Pparam(uesp_taper_value, uesp_taper_range, uesp_taper_dstep, uesp_taper_vary, uesp_taper_defined); + uesp_temp = Pparam(uesp_temp_value, uesp_temp_range, uesp_temp_dstep, uesp_temp_vary, uesp_temp_defined); + } + + this->delta_phase = delta_phase; + this->nlat1f = nlat1f, this->nlat2f = nlat2f, this->nlat1c = nlat1c, this->nlat2c = nlat2c; + this->npole = npole; this->nlatfill = nlatfill, this->nlngfill = nlngfill; + this->lfudge = lfudge, this->llo = llo, this->lhi = lhi; + this->phase1 = phase1, this->phase2 = phase2, this->nrad = nrad; + this->wavelength = wavelength; + this->roche1 = roche1, this->roche2 = roche2, this->eclipse1 = eclipse1, this->eclipse2 = eclipse2; + this->glens1 = glens1, this->use_radii = use_radii; + this->tperiod = tperiod, this->gdark_bolom1 = gdark_bolom1, this->gdark_bolom2 = gdark_bolom2; + this->mucrit1 = mucrit1, this->mucrit2 = mucrit2; + this->mirror = mirror, this->add_disc = add_disc, this->opaque = opaque, this->add_spot = add_spot, this->nspot = nspot; + this->iscale = iscale; + if(glens1 && roche1) + throw Lcurve_Error("For reasons of simplicity, glens1 = 1 and roche1 = 1 are not simultaneously allowed"); + std::string slimb1 = pslimb1; + std::string slimb2 = pslimb2; + //if(strcmp(pslimb1, "Poly") == 0){ + if(slimb1 == "Poly"){ + limb1 = LDC::POLY; + }else if(slimb1 == "Claret"){ + limb1 = LDC::CLARET; + }else{ + throw Lcurve_Error("Could not recognize the value of slimb1 = " + + slimb1 + "; should be 'Poly' or 'Claret'"); + } + if(slimb2 == "Poly"){ + limb2 = LDC::POLY; + }else if(slimb2 == "Claret"){ + limb2 = LDC::CLARET; + }else{ + throw Lcurve_Error("Could not recognize the value of slimb2 = " + + slimb2 + "; should be 'Poly' or 'Claret'"); + } + +} + + + /** Works out the number of variable physical parameters */ int Lcurve::Model::nvary() const { int n = 0; From aa1bdefc9b13b6db5eabc6bbbc989ee0213cf0f7 Mon Sep 17 00:00:00 2001 From: lidihei <912219465@qq.com> Date: Wed, 19 Oct 2022 14:41:51 +0800 Subject: [PATCH 02/12] add libpylcurve.cc --- src/libpylcurve.cc | 606 ++++++++++++++++++++++++++++++++++++++ src/pylight_curve_comp.cc | 312 ++++++++++++++++++++ 2 files changed, 918 insertions(+) create mode 100644 src/libpylcurve.cc create mode 100644 src/pylight_curve_comp.cc diff --git a/src/libpylcurve.cc b/src/libpylcurve.cc new file mode 100644 index 0000000..c825d05 --- /dev/null +++ b/src/libpylcurve.cc @@ -0,0 +1,606 @@ +/* + +!!begin +!!title Computes the light-curve of a sphere and a Roche distorted star +!!author Jiao Li +!!created 08 Oct 2022 +!!descr Computes the light-curve of a sphere and a Roche distorted star +!!index lroche +!!lib libpylcurve +!!css style.css +!!class Model +!!head1 Computes the light-curve of a sphere and a Roche distorted star + +!!emph{lroche} computes the light curve equivalent to a model of a sphere and +a Roche-distorted star to model a white dwarf/main-sequence binary such as NN +Ser and can optionally include a disc and bright-spot as well. It includes +reprocessing and eclipses. The reprocessing is computed using the simple +addition of fluxes method, assuming that the irradiating star can be treated +as a point, and not including any 'back heating'. Phase 0 is defined as the +point when star 1 is furthest from the observer. The star grids can de adapted +somewhat to be focussed on eclipses and long exposures trapezoidally sub-divided +to model smearing. Other physics included in !!emph{lroche}: Doppler beaming, +gravitational lensing, Roemer time delays, asynchronous rotation of the stellar components. +Physics not included: eccentric orbits! + +!!emph{lroche} is also designed as a test routine and so can add gaussian +noise to the result. It can either read a data file and compute at exactly the +right times or it can compute the light curve over a regularly spaced +sequence. This file is the main documentation for the data and model file +formats. + +!!head2 Invocation + +lroche model data [time1 time2 ntime expose] noise seed nfile [output] (device) [(roff) scale [ssfac]/[sstar1 sstar2 sdisc spot]] + +!!head2 Arguments + +!!table +!!arg{model}{File of parameters specifying the parameters. See below for a full description of these.} +!!arg{data}{Data file as a template to compute the fit. 'none' to indicate that +you do not have one and will define times instead. In the case of real data, the fit +will be scaled automatically to minimise chi**2. Data files requires six columns which +are: mid-exposure times, exposure times (same units as times), fluxes, uncertainties on fluxes, weight factors and +finally integers to represent the number of subdivisions to model finite exposure effects. The weight factors +allow you to change the weighting of a point to the overall chi**2 without changing its uncertainty.} +!!arg{time1}{If data = 'none', first time to compute.} +!!arg{time2}{If data = 'none', last time to compute.} +!!arg{ntime}{If data = 'none', number of times to compute.} +!!arg{expose}{If data = 'none', length of exposure} +!!arg{noise}{If data = 'none', amount of noise to add to the results, RMS. In the case of real data, +this is used a multiplier of the real error bars.} +!!arg{seed}{Seed integer} +!!arg{nfile}{Number of files to store, each with its own noise; only the first will be plotted if +nfile > 1. 0 is possible in which case there will be no output.} +!!arg{output}{File to save the results in the form of rows each with +time, exposure time, flux and uncertainty. If nfile > 1, then the files +will have a number added, as in data001} +!!arg{device}{Plot device to use; 'none' or 'null' to ignore} +!!arg{roff}{Offset to add to the residuals when plotting a fit to data} +!!arg{scale}{Autoscale or not. If true, then the scale factors will be determined by minimisation of chi**2 using +linear scaling factors.} +!!arg{ssfac}{If scale=false, and a single global scaling factor is set, this is its value. Depending on the units of the fluxes, ssfac has the following interpretation. If ssfac is set = +(a/d)^2 (a = semi-major axis, d = distance), the result returned will be the flux density at Earth (fnu) in SI units (W/m^2/Hz). A milliJansky = 10^-29 W/m^2/Hz, so if you are fitting mJy +fluxes, the value of ssfac returned = 10^29 (a/d)^2. If you express a in solar radii and d in parsecs, and fit mJy fluxes, then ssfac = 5.0875 x 10^13 (a/d)^2.} +!!arg{sstar1}{If scale=false, scale factor for star 1} +!!arg{sstar2}{If scale=false, scale factor for star 2} +!!arg{sdisc}{If scale=false, scale factor for disc} +!!arg{sspot}{If scale=false, scale factor for spot} + +!!table + +!!head2 Parameter file + +The model parameters come in two types, physical and computational. 'Physical' in this case are ones which have an initial +value, a range of plausible variation and a step size for derivative computation. 'Computational' are ones which cannot +be varied and just have a value. !!emph{lroche} ignores the variable/non-variable distinction which are for the fitting +routines such as !!ref{levmarq.html}{levmarq}. The parameter file consists of a series of lines such as: +
+q      =  0.12  0.01 0.0001 0 1
+iangle =  82    2    0.01   1 1
+r1     =  0.17  0.05 0.001  1 1
+.
+.
+.
+etc
+
+delta_phase = 1.e-7 .  .  .  etc
etc. The above would mean that q is +not to be varied (0 as the penultimate value), but iangle and r1 are. For +!!ref{simplex.html}{simplex}, !!ref{genetic.html}{genetic}, +!!ref{simann.html}{simann} and !!ref{powell.html}{powell} the second parameter +specifies the range over which to vary the respective parameter (i.e. 0.01, 2 +and 0.05 in this case), while the third parameter is used by +!!ref{levmarq.html}{levmarq} to compute numerical derivatives using finite +differences. Be careful to set this small enough to give accurate derivatives +but not so small that roundoff will be problematic. It should, at a minimum, be +smaller than the uncertainty in any parameter. The final integer is a relatively recent addition which specifies whether a given parameter is used at all and is designed for the star spot parameters. + +Here is the full list of parameters. + +!!head2 Physical parameters + +!!head3 Binary and stars + +!!table +!!arg{q}{Mass ratio, q = M2/M1} +!!arg{iangle}{Inclination angle, degrees} +!!arg{r1}{Radius of star 1, scaled by the binary separation} +!!arg{r2}{Radius of star 2, scaled by the binary separation. The radius is measured along the line of centres towards star 1. Set = -1 and hold fixed for +Roche lobe filling stars.} +!!arg{cphi3}{Third contact phase (star 1 starting to emerge from eclipse). This is an alternative way to specify the radii, based on a +spherical approximation fot the two stars, i.e. unless the stars are spherical, it is not quite the true third contact. The radii will +be computed from the contact phases according to the two equations r2+r1 = sqrt(1 - sin^2 i cos^2 (2*pi*cphi4)) and +r2-r1 = sqrt(1 - sin^2 i cos^2 (2*pi*cphi3)). The radii returned are precise, just the interpretation as contact phases that is not +precise. cphi3 and cphi4 need the boolean use_radii set to 0 to enabled. The reason for using them is to help with MCMC iterations as +they prevent the nasty curved correlation between r1, r2 and i. This can save a huge amount of CPU time.} +!!arg{cphi4}{Fourth contact phase, star 1 fully emerged from eclipse. See cphi3 for details.} +!!arg{spin1}{This is the ratio of the spin frequency of star 1 to the orbital frequency. In this case a modified form of the Roche potential is used for star 1} +!!arg{spin2}{This is the ratio of the spin frequency of star 2 to the orbital frequency. In this case a modified form of the Roche potential is used for star 2} +!!arg{t1}{Temperature of star 1, Kelvin. This is really a substitute for surface brightness which is set assuming a black-body given this parameter. If it was not for irradiation that would be exactly what this is, a one-to-one replacement for surface brightness. Irradiation however introduces bolometric luminosities effectively and breaks the direct link. Some would then argue that one must use model atmospheres except at the moment irradiated model atmosphere are in their infancy.} +!!arg{t2}{Temperature of star 2, Kelvin. Set < 0 in order that it does not get scaled when using the iscale parameter.} +!!arg{ldc1_1, etc}{Limb darkening for stars is quite hard to specify precisely. Here we adopt a 4 coefficient +approach which can either represent a straighforward polynimal expanion of the form I(mu) = 1 - \sum_i a_i (1-mu)^i, +or rather better in some cases Claret's 4-coefficient formula I(mu) = 1 - \sum_i a_i (1 - mu^(i/2)) (i=1 to 4). +You specify these by supplying the 4 coefficients for each star (which for form's sake are potentially variable but +you would probably be unwise to let them be free) and later on a parameter to say whether it is the polynomial or +Claret's representation. The polynomial allows one to use linear and quadratic limb darkening amongst others by setting +the upper coefficients = 0. ldc1_1 is the first coefficient of star 1, ldc1_2 is the second, etc, while ldc2_1 +is the first coefficient for star 2 etc. See limb1, limb2, mucrit1, mucrit2 below.} +!!arg{velocity_scale}{Velocity scale, sum of unprojected orbital speeds, used for accounting for +Doppler beaming and gravitational lensing. On its own this makes little difference to the light curve, so you should +not usually let it be free, but you might want to if you have independent K1 or K2 information which you can apply as part of a prior.} +!!arg{beam_factor1}{The factor to use for Doppler beaming from star 1. This corresponds to the factor (3-alpha) +that multiplies -v_r/c in the standard beaming formula where alpha is related to the spectral shape. Use of this parameter +requires the velocity_scale to be set.} +!!arg{beam_factor2}{The factor to use for Doppler beaming from star 2. This corresponds to the factor (3-alpha) +that multiplies -v_r/c in the standard beaming formula where alpha is related to the spectral shape. Use of this parameter +requires the velocity_scale to be set.} +!!table + +!!head3 General + +!!table +!!arg{t0}{Zero point of ephemeris, marking time of mid-eclipse (or in general +superior conjunction) of star 1, same units as times.} +!!arg{period}{Orbital period, same units as times.} +!!arg{pdot}{Quadratic coefficient of ephemeris, same units as times} +!!arg{deltat}{Time shift between the primary and secondary eclipses to allow for small eccentricities and Roemer delays in the orbit. The +sign is defined such that deltat > 0 implies that the secondary eclipse suffers a delay compared to the primary compared to precisely 0.5 +difference. deltat < 0 implies the secondary eclipse comes a little earlier than expected. Assuming that the "primary eclipse" is the eclipse +of star 1, then, using the same sign convention, the Roemer delay is given by = P*(K1-K2)/(Pi*c) where P is the orbital period, K1 and K2 are +the usual projected radial velocity semi-amplitudes Pi = 3.14159.., and c = speed of light. See Kaplan (2010) for more details. The delay +is implemented by adjusting the orbital phase according to phi' = phi + (deltat/2/P)*(cos(2*Pi*phi)-1), i.e. there is no change at primary eclipse +but a delay of -deltat/P by the secondary eclipse.} +!!arg{gravity_dark}{Gravity darkening coefficient. Only matters for the Roche distorted case, but is prompted for always. +There are two alternatives for this. In the standard old method, the temperatures on the stars are set equal to t2*(g/gr)**gdark +where g is the gravity at a given point and gr is the gravity at the point furthest from the primary (the 'backside' of the secondary). +For a convectuive atmosphere, 0.08 is the usual value while 0.25 is the number for a radiative atmosphere. This is translated into +intensity using a blackbody approx. If you want to bypass the BB approx and invoke a direct relation flux ~ (g/gr)**gdark relation +you should set gdark_bolom (see below) to 0 (false.)} +!!arg{absorb}{The fraction of the irradiating flux from star 1 absorbed by star 2} +!!arg{slope, quad, cube}{Fudge factors to help cope with any trends in the data +as a result of e.g. airmass effects. The fit is multiplied by +(1+x*(slope+x*(quad+x*cube))) where x is the time scaled so that it varies +from -1 to 1 from start to end of the data. One should expect these number +to have absolute value << 1.} +!!arg{third}{Third light contribution. Simply adds to whatever flux is +calculated and will be subject to auto-scaling like other flux. It only +applies if global scaling rather than individual component scaling is used. +Third light is assumed strictly constant and is not subject to the slope, +quad, cube parameters.} +!!table + +!!head3 Spots + +One spot allowed on each star (with some expectation that the +number may be increased if need be in the future): + +!!table +!!arg{stsp11_long}{Longitude (degrees) of spot 1 on star 1, relative to +meridian defined by line of centres} +!!arg{stsp11_lat}{Latitude (degrees) of spot 1 on star 1} +!!arg{stsp11_lat}{FWHM (degrees) of spot 1 on star 1, as seen from its centre +of mass. Spot has gaussian distribution of temperature.} +!!arg{stsp11_tcen}{Central temp (K) of spot 1 on star 1} +!!arg{stsp21_long}{Longitude (degrees) of spot 1 on star 2, relative to +meridian defined by line of centres} +!!arg{stsp21_lat}{Latitude (degrees) of spot 1 on star 2} +!!arg{stsp21_lat}{FWHM (degrees) of spot 1 on star 2, as seen from its centre +of mass. Spot has gaussian distribution of temperature.} +!!arg{stsp21_tcen}{Central temp (K) of spot 1 on star 2} +!!table + +!!head3 Disc + +!!table +!!arg{rdisc1}{Inner radius of azimuthally symmetric disc. Set = -1 to set it equal to r1 (it should not be allowed to vary in this case)} +!!arg{rdisc2}{Outer radius of azimuthally symmetric disc. Set = -1 and hold fixed to clamp this to equal the bright spot radius.} +!!arg{height_disc}{Half height of disc at radius = 1. The height varies as a power law of radius} +!!arg{beta_disc}{Exponent of power law in radius of disc. Should be >= 1 to make concave disc; convex will not eclipse +properly.} +!!arg{temp_disc}{Temperature of outer part of disc. This is little more than a flux normalisation parameter but +it is easier to think in terms of temperature} +!!arg{texp_disc}{Exponent of surface brightness (NB: not temperature) over disc} +!!arg{lin_limb_disc}{Linear limb darkening coefficient of the disc} +!!arg{quad_limb_disc}{Quadratic limb darkening coefficient of the disc} +!!arg{temp_edge}{Temperature at perpendicular edge of disc. Irradiation from the secondary is allowed so you should think of a bright rim at primary eclipse. Limb darkeining parameters of the +disc are applied} +!!arg{absorb_edge}{Amount of secondary flux absorbed and reprocessed. This effect should lead +to a sinusoidal variation with flux maximum at orbital phase 0.5. It was introduced to model +a possible accreting sdO/WD system discovered by Thomas Kupfer} + +!!table + +!!head3 Bright-spot + +!!table +!!arg{radius_spot}{Distance from accretor of bright-spot (units of binary separation).} +!!arg{length_spot}{Length scale of spot (units of binary separation).} +!!arg{height_spot}{Height of spot (units of binary separation). This is only a normalisation constant.} +!!arg{expon_spot}{Spot is modeled as x**n*exp(-(x/l)**m). This parameter specifies the exponent 'n'} +!!arg{epow_spot}{This is the exponent m in the above expression} +!!arg{angle_spot}{This is the angle made by the line of elements of the spot measured in the direction of binary motion relative to +the rim of the disc so that the "standard" value should be 0.} +!!arg{yaw_spot}{Allows the spot elements effectively to beam their light away from the perpendicular to the line of elements. +Measured as an angle in the same sense as angle_spot. 0 means standard perpendicular beaming.} +!!arg{temp_spot}{Normalises the surface brightness of the spot.} +!!arg{tilt_spot}{Allows spot to be other than perpendicular to the disc. 90 = perpendicular. If less than 90 then the +spot is visible for more than half a cycle.} +!!arg{cfrac_spot}{The fraction of the spot taken to be equally visible at all phases, i.e. pointing upwards.} +!!table + +!!head2 Computational parameters +!!table +!!arg{delta_phase}{Accuracy in phase of eclipse computations. This determines the +accuracy of any Roche computations. Example: 1.e-7} +!!arg{nlat1f}{The number of latitudes for star 1's fine grid. This is used around the phase of primary eclipse (i.e. the +eclipse of star 1} +!!arg{nlat1c}{The number of latitudes for star 1's coarse grid. This is used away from primary eclipse.} +!!arg{nlat2f}{The number of latitudes for star 2's fine grid. This is used around the phase of secondary eclipse.} +!!arg{nlat2c}{The number of latitudes for star 2's coarse grid. This is used away from secondary eclipse.} +!!arg{npole}{True to set North pole of grid to the genuine stellar NP rather than substellar points. This is probably a good idea +when modelling well detached binaries, especially with extreme radius ratios because then it allows one to concentrate points +over a band of latitudes using the next two parameters} +!!arg{nlatfill}{Extra number of points to insert per normal latitude strip along the path of star 1 as it transits star 2. This is +designed to help tough extreme radius ratio cases. Take care to look at the resulting grid with visualise as the exact latitude range +chosen is a little approximate. This is only enabled if npole since only then do the latitude strips more-or-less line up with the movement +of the star.} +!!arg{nlngfill}{Extra number of points to insert per normal longitude strip along the path of star 1 as it transits star 2. This is +designed to help tough extreme radius ratio cases. Take care to look at the resulting grid with visualise as the exact latitude range +chosen is a little approximate.} + +!!arg{lfudge}{The fine-grid latitude strip is computed assuming both stars are spherical. To allow for departures from this, this parameter +allows one to increase the latitude limits both up and down by an amount specified in degrees. Use the program !!ref{visualise.html}{visualise} +to judge how large this should be. However, one typically would like to avoid lfudge > 30*r1/r2 as that could more than double the width of the strip.} + +!!arg{phase1}{this defines when star 1's fine grid is used abs(phase) < phase1. Thus phase1 = 0.05 will restrict +the fine grid use to phase 0.95 to 0.05.} +!!arg{llo, lhi}{These are experimental. They allow the user to fix the latitude limits of the fine strip which might be useful in preventing chi**2 variations +caused by variable grids. The values need to reflect the likely range of inclinations and can only really be set by trial and error using visualise. They are in degrees +following the usual convention for latitude on Earth. Set llo high and lhi low to stop them having any effect.} +!!arg{phase1}{this defines when star 1's fine grid is used abs(phase) < phase1. Thus phase1 = 0.05 will restrict +the fine grid use to phase 0.95 to 0.05.} +!!arg{phase2}{this defines when star 2's fine grid is used phase2 until 1-phase2. Thus phase2 = 0.45 will restrict +the fine grid use to phase 0.55 to 0.55.} +!!arg{nrad}{The number of radial strips over the disc} +!!arg{wavelength}{Wavelength (nm)} +!!arg{roche1}{Account for Roche distortion of star 1 or not} +!!arg{roche2}{Account for Roche distortion of star 2 or not} +!!arg{eclipse1}{Account for the eclipse of star 1 or not} +!!arg{eclipse2}{Account for the eclipse of star 2 or not} +!!arg{glens1}{Account for gravitational lensing by star 1. If you use this roche1 must be = 0 and the velocity_scale} +!!arg{use_radii}{If set = 1, the parameters r1 and r2 will be used to set the radii directly. If not, the third and fourth contact phases, +cphi3 and cphi4, will be used instead (see description for cphi3 for details).} +!!arg{tperiod}{The true orbital period in days. This is required, along with velocity_scale, if gravitational lensing is +being applied to calculate proper dimensions in the system.} +!!arg{gdark_bolom}{True if the gravity darkening coefficient represents the bolometric value where T is proportional to +gravity to the power set by the coefficient. This is translated to flux variations using the black-body approximation. +If False, it represents a filter-integrated value 'y' coefficient such that the flux depends upon the gravity to the power 'y'. +This is itself an approximation and ideally should replaced by a proper function of gravity, but is probably good enough for +most purposes. Please see gravity_dark.} +!!arg{mucrit1}{Critical value of mu on star 1 below which intensity is assumed to be zero. This is to allow one to represent +Claret and Hauschildt's (2004) results where I(mu) drops steeply for mu < 0.08 or so. WARNING: this option is dangerous. I would normally +advise setting it = 0 unless you really know what you are doing as it leads to discontinuities.} +!!arg{mucrit2}{Critical value of mu on star 2 below which intensity is assumed to be zero. See comments on mucrit1 for more.} +!!arg{limb1}{String, either 'Poly' or 'Claret' determining the type of limb darkening law. See comments on ldc1_1 above.} +!!arg{limb2}{String, either 'Poly' or 'Claret' determining the type of limb darkening law. See comments on ldc1_1 above.} +!!arg{mirror}{Add any light not reprocessed in as if star reflected it or not as a crude approximation to the +effet of gray scattering} +!!arg{add_disc}{Add a disc or not} +!!arg{opaque}{Make disc opaque or not} +!!arg{iscale}{Individually scale the separate components or not. If set the +each component, star 1, star 2, disc and bright spot will be individually +scaled to minimise chi**2. Otherwise a single overall factor will be computed. +NB If you set this parameter then all temperature parameters (white dwarf, +secondary, disc and bright spot) must be held fixed otherwise near-total +degeneracy will result. The only reason it is not total is because of +reflection effect from irradiation of the secondary by the white dwarf, but +this is often very feeble and will not help, so, you have been warned. +Scaling should in general lead to faster convergence than not scaling. +You may have some cases where you do not want to include any secondary +star component. You can do this by setting t2 < 0. Note that if this is set +true, then the third light parameter will be ignored.} +!!table + +!!end +* This file is modified by lijiao +*/ + +#include +#include +#include +#include +#include "cpgplot.h" +#include "trm/subs.h" +#include "trm/plot.h" +#include "trm/vec3.h" +#include "trm/input.h" +#include "trm/format.h" +#include "trm/roche.h" +#include "trm/lcurve.h" + +int Lcurve::Fobj::neval = 0; +double Lcurve::Fobj::chisq_min; +Subs::Buffer1D Lcurve::Fobj::scale_min; + +// Main program +extern "C"{ + void pylcurve(double *times, double *exposes, int *ndivs, int Tsize, + double *calc, double *lcstar1, double *lcdisc, + double *lcedge, double *lcspot, double *lcstar2, double *wdwarflogrv, + //Binary and stars + double q_value, double q_range, double q_dstep, bool q_vary, bool q_defined, + double iangle_value, double iangle_range, double iangle_dstep, bool iangle_vary, bool iangle_defined, + double r1_value, double r1_range, double r1_dstep, bool r1_vary, bool r1_defined, + double r2_value, double r2_range, double r2_dstep, bool r2_vary, bool r2_defined, + double cphi3_value, double cphi3_range, double cphi3_dstep, bool cphi3_vary, bool cphi3_defined, + double cphi4_value, double cphi4_range, double cphi4_dstep, bool cphi4_vary, bool cphi4_defined, + double spin1_value, double spin1_range, double spin1_dstep, bool spin1_vary, bool spin1_defined, + double spin2_value, double spin2_range, double spin2_dstep, bool spin2_vary, bool spin2_defined, + double t1_value, double t1_range, double t1_dstep, bool t1_vary, bool t1_defined, + double t2_value, double t2_range, double t2_dstep, bool t2_vary, bool t2_defined, + double ldc1_1_value, double ldc1_1_range, double ldc1_1_dstep, bool ldc1_1_vary, bool ldc1_1_defined, + double ldc1_2_value, double ldc1_2_range, double ldc1_2_dstep, bool ldc1_2_vary, bool ldc1_2_defined, + double ldc1_3_value, double ldc1_3_range, double ldc1_3_dstep, bool ldc1_3_vary, bool ldc1_3_defined, + double ldc1_4_value, double ldc1_4_range, double ldc1_4_dstep, bool ldc1_4_vary, bool ldc1_4_defined, + double ldc2_1_value, double ldc2_1_range, double ldc2_1_dstep, bool ldc2_1_vary, bool ldc2_1_defined, + double ldc2_2_value, double ldc2_2_range, double ldc2_2_dstep, bool ldc2_2_vary, bool ldc2_2_defined, + double ldc2_3_value, double ldc2_3_range, double ldc2_3_dstep, bool ldc2_3_vary, bool ldc2_3_defined, + double ldc2_4_value, double ldc2_4_range, double ldc2_4_dstep, bool ldc2_4_vary, bool ldc2_4_defined, + double velocity_scale_value, double velocity_scale_range, double velocity_scale_dstep, bool velocity_scale_vary, bool velocity_scale_defined, + double beam_factor1_value, double beam_factor1_range, double beam_factor1_dstep, bool beam_factor1_vary, bool beam_factor1_defined, + double beam_factor2_value, double beam_factor2_range, double beam_factor2_dstep, bool beam_factor2_vary, bool beam_factor2_defined, + //General + double t0_value, double t0_range, double t0_dstep, bool t0_vary, bool t0_defined, + double period_value, double period_range, double period_dstep, bool period_vary, bool period_defined, + double pdot_value, double pdot_range, double pdot_dstep, bool pdot_vary, bool pdot_defined, + double deltat_value, double deltat_range, double deltat_dstep, bool deltat_vary, bool deltat_defined, + double gravity_dark1_value, double gravity_dark1_range, double gravity_dark1_dstep, bool gravity_dark1_vary, bool gravity_dark1_defined, + double gravity_dark2_value, double gravity_dark2_range, double gravity_dark2_dstep, bool gravity_dark2_vary, bool gravity_dark2_defined, + double absorb_value, double absorb_range, double absorb_dstep, bool absorb_vary, bool absorb_defined, + double slope_value, double slope_range, double slope_dstep, bool slope_vary, bool slope_defined, + double quad_value, double quad_range, double quad_dstep, bool quad_vary, bool quad_defined, + double cube_value, double cube_range, double cube_dstep, bool cube_vary, bool cube_defined, + double third_value, double third_range, double third_dstep, bool third_vary, bool third_defined, + // Star spots + double stsp11_long_value, double stsp11_long_range, double stsp11_long_dstep, bool stsp11_long_vary, bool stsp11_long_defined, + double stsp11_lat_value, double stsp11_lat_range, double stsp11_lat_dstep, bool stsp11_lat_vary, bool stsp11_lat_defined, + double stsp11_fwhm_value, double stsp11_fwhm_range, double stsp11_fwhm_dstep, bool stsp11_fwhm_vary, bool stsp11_fwhm_defined, + double stsp11_tcen_value, double stsp11_tcen_range, double stsp11_tcen_dstep, bool stsp11_tcen_vary, bool stsp11_tcen_defined, + double stsp12_long_value, double stsp12_long_range, double stsp12_long_dstep, bool stsp12_long_vary, bool stsp12_long_defined, + double stsp12_lat_value, double stsp12_lat_range, double stsp12_lat_dstep, bool stsp12_lat_vary, bool stsp12_lat_defined, + double stsp12_fwhm_value, double stsp12_fwhm_range, double stsp12_fwhm_dstep, bool stsp12_fwhm_vary, bool stsp12_fwhm_defined, + double stsp12_tcen_value, double stsp12_tcen_range, double stsp12_tcen_dstep, bool stsp12_tcen_vary, bool stsp12_tcen_defined, + double stsp13_long_value, double stsp13_long_range, double stsp13_long_dstep, bool stsp13_long_vary, bool stsp13_long_defined, + double stsp13_lat_value, double stsp13_lat_range, double stsp13_lat_dstep, bool stsp13_lat_vary, bool stsp13_lat_defined, + double stsp13_fwhm_value, double stsp13_fwhm_range, double stsp13_fwhm_dstep, bool stsp13_fwhm_vary, bool stsp13_fwhm_defined, + double stsp13_tcen_value, double stsp13_tcen_range, double stsp13_tcen_dstep, bool stsp13_tcen_vary, bool stsp13_tcen_defined, + double stsp21_long_value, double stsp21_long_range, double stsp21_long_dstep, bool stsp21_long_vary, bool stsp21_long_defined, + double stsp21_lat_value, double stsp21_lat_range, double stsp21_lat_dstep, bool stsp21_lat_vary, bool stsp21_lat_defined, + double stsp21_fwhm_value, double stsp21_fwhm_range, double stsp21_fwhm_dstep, bool stsp21_fwhm_vary, bool stsp21_fwhm_defined, + double stsp21_tcen_value, double stsp21_tcen_range, double stsp21_tcen_dstep, bool stsp21_tcen_vary, bool stsp21_tcen_defined, + double stsp21_long_value, double stsp21_long_range, double stsp21_long_dstep, bool stsp21_long_vary, bool stsp21_long_defined, + double stsp21_lat_value, double stsp21_lat_range, double stsp21_lat_dstep, bool stsp21_lat_vary, bool stsp21_lat_defined, + double stsp21_fwhm_value, double stsp21_fwhm_range, double stsp21_fwhm_dstep, bool stsp21_fwhm_vary, bool stsp21_fwhm_defined, + double stsp21_tcen_value, double stsp21_tcen_range, double stsp21_tcen_dstep, bool stsp21_tcen_vary, bool stsp21_tcen_defined, + double uesp_long1_value, double uesp_long1_range, double uesp_long1_dstep, bool uesp_long1_vary, bool uesp_long1_defined, + double uesp_long2_value, double uesp_long2_range, double uesp_long2_dstep, bool uesp_long2_vary, bool uesp_long2_defined, + double uesp_lathw_value, double uesp_lathw_range, double uesp_lathw_dstep, bool uesp_lathw_vary, bool uesp_lathw_defined, + double uesp_taper_value, double uesp_taper_range, double uesp_taper_dstep, bool uesp_taper_vary, bool uesp_taper_defined, + double uesp_temp_value, double uesp_temp_range, double uesp_temp_dstep, bool uesp_temp_vary, bool uesp_temp_defined, + // disc + double rdisc1_value, double rdisc1_range, double rdisc1_dstep, bool rdisc1_vary, bool rdisc1_defined, + double rdisc2_value, double rdisc2_range, double rdisc2_dstep, bool rdisc2_vary, bool rdisc2_defined, + double height_disc_value, double height_disc_range, double height_disc_dstep, bool height_disc_vary, bool height_disc_defined, + double beta_disc_value, double beta_disc_range, double beta_disc_dstep, bool beta_disc_vary, bool beta_disc_defined, + double temp_disc_value, double temp_disc_range, double temp_disc_dstep, bool temp_disc_vary, bool temp_disc_defined, + double texp_disc_value, double texp_disc_range, double texp_disc_dstep, bool texp_disc_vary, bool texp_disc_defined, + double lin_limb_disc_value, double lin_limb_disc_range, double lin_limb_disc_dstep, bool lin_limb_disc_vary, bool lin_limb_disc_defined, + double quad_limb_disc_value, double quad_limb_disc_range, double quad_limb_disc_dstep, bool quad_limb_disc_vary, bool quad_limb_disc_defined, + double temp_edge_value, double temp_edge_range, double temp_edge_dstep, bool temp_edge_vary, bool temp_edge_defined, + double absorb_edge_value, double absorb_edge_range, double absorb_edge_dstep, bool absorb_edge_vary, bool absorb_edge_defined, + //Bright-spot + double radius_spot_value, double radius_spot_range, double radius_spot_dstep, bool radius_spot_vary, bool radius_spot_defined, + double length_spot_value, double length_spot_range, double length_spot_dstep, bool length_spot_vary, bool length_spot_defined, + double height_spot_value, double height_spot_range, double height_spot_dstep, bool height_spot_vary, bool height_spot_defined, + double expon_spot_value, double expon_spot_range, double expon_spot_dstep, bool expon_spot_vary, bool expon_spot_defined, + double epow_spot_value, double epow_spot_range, double epow_spot_dstep, bool epow_spot_vary, bool epow_spot_defined, + double angle_spot_value, double angle_spot_range, double angle_spot_dstep, bool angle_spot_vary, bool angle_spot_defined, + double yaw_spot_value, double yaw_spot_range, double yaw_spot_dstep, bool yaw_spot_vary, bool yaw_spot_defined, + double temp_spot_value, double temp_spot_range, double temp_spot_dstep, bool temp_spot_vary, bool temp_spot_defined, + double tilt_spot_value, double tilt_spot_range, double tilt_spot_dstep, bool tilt_spot_vary, bool tilt_spot_defined, + double cfrac_spot_value, double cfrac_spot_range, double cfrac_spot_dstep, bool cfrac_spot_vary, bool cfrac_spot_defined, + // Computational parameters + double delta_phase, int nlat1f, int nlat2f, int nlat1c, int nlat2c, bool npole, + int nlatfill, int nlngfill, double lfudge, double llo, double lhi, double phase1, double phase2, int nrad, double wavelength, + bool roche1, bool roche2, bool eclipse1, bool eclipse2, bool glens1, bool use_radii, + double tperiod, double gdark_bolom1, double gdark_bolom2, double mucrit1, double mucrit2, + const char* pslimb1, const char* pslimb2, bool mirror, bool add_disc, bool opaque, bool add_spot, int nspot, bool iscale, bool info, + int parallel_threshold + ){ + try{ + //added by lijiao + Lcurve::Model model( // Binary and stars + q_value, q_range, q_dstep, q_vary, q_defined, + iangle_value, iangle_range, iangle_dstep, iangle_vary, iangle_defined, + r1_value, r1_range, r1_dstep, r1_vary, r1_defined, + r2_value, r2_range, r2_dstep, r2_vary, r2_defined, + cphi3_value, cphi3_range, cphi3_dstep, cphi3_vary, cphi3_defined, + cphi4_value, cphi4_range, cphi4_dstep, cphi4_vary, cphi4_defined, + spin1_value, spin1_range, spin1_dstep, spin1_vary, spin1_defined, + spin2_value, spin2_range, spin2_dstep, spin2_vary, spin2_defined, + t1_value, t1_range, t1_dstep, t1_vary, t1_defined, + t2_value, t2_range, t2_dstep, t2_vary, t2_defined, + ldc1_1_value, ldc1_1_range, ldc1_1_dstep, ldc1_1_vary, ldc1_1_defined, + ldc1_2_value, ldc1_2_range, ldc1_2_dstep, ldc1_2_vary, ldc1_2_defined, + ldc1_3_value, ldc1_3_range, ldc1_3_dstep, ldc1_3_vary, ldc1_3_defined, + ldc1_4_value, ldc1_4_range, ldc1_4_dstep, ldc1_4_vary, ldc1_4_defined, + ldc2_1_value, ldc2_1_range, ldc2_1_dstep, ldc2_1_vary, ldc2_1_defined, + ldc2_2_value, ldc2_2_range, ldc2_2_dstep, ldc2_2_vary, ldc2_2_defined, + ldc2_3_value, ldc2_3_range, ldc2_3_dstep, ldc2_3_vary, ldc2_3_defined, + ldc2_4_value, ldc2_4_range, ldc2_4_dstep, ldc2_4_vary, ldc2_4_defined, + velocity_scale_value, velocity_scale_range, velocity_scale_dstep, velocity_scale_vary, velocity_scale_defined, + beam_factor1_value, beam_factor1_range, beam_factor1_dstep, beam_factor1_vary, beam_factor1_defined, + beam_factor2_value, beam_factor2_range, beam_factor2_dstep, beam_factor2_vary, beam_factor2_defined, + //General + t0_value, t0_range, t0_dstep, t0_vary, t0_defined, + period_value, period_range, period_dstep, period_vary, period_defined, + pdot_value, pdot_range, pdot_dstep, pdot_vary, pdot_defined, + deltat_value, deltat_range, deltat_dstep, deltat_vary, deltat_defined, + gravity_dark1_value, gravity_dark1_range, gravity_dark1_dstep, gravity_dark1_vary, gravity_dark1_defined, + gravity_dark2_value, gravity_dark2_range, gravity_dark2_dstep, gravity_dark2_vary, gravity_dark2_defined, + absorb_value, absorb_range, absorb_dstep, absorb_vary, absorb_defined, + slope_value, slope_range, slope_dstep, slope_vary, slope_defined, + quad_value, quad_range, quad_dstep, quad_vary, quad_defined, + cube_value, cube_range, cube_dstep, cube_vary, cube_defined, + third_value, third_range, third_dstep, third_vary, third_defined, + // Star spots + stsp11_long_value, stsp11_long_range, stsp11_long_dstep, stsp11_long_vary, stsp11_long_defined, + stsp11_lat_value, stsp11_lat_range, stsp11_lat_dstep, stsp11_lat_vary, stsp11_lat_defined, + stsp11_fwhm_value, stsp11_fwhm_range, stsp11_fwhm_dstep, stsp11_fwhm_vary, stsp11_fwhm_defined, + stsp11_tcen_value, stsp11_tcen_range, stsp11_tcen_dstep, stsp11_tcen_vary, stsp11_tcen_defined, + stsp12_long_value, stsp12_long_range, stsp12_long_dstep, stsp12_long_vary, stsp12_long_defined, + stsp12_lat_value, stsp12_lat_range, stsp12_lat_dstep, stsp12_lat_vary, stsp12_lat_defined, + stsp12_fwhm_value, stsp12_fwhm_range, stsp12_fwhm_dstep, stsp12_fwhm_vary, stsp12_fwhm_defined, + stsp12_tcen_value, stsp12_tcen_range, stsp12_tcen_dstep, stsp12_tcen_vary, stsp12_tcen_defined, + stsp13_long_value, stsp13_long_range, stsp13_long_dstep, stsp13_long_vary, stsp13_long_defined, + stsp13_lat_value, stsp13_lat_range, stsp13_lat_dstep, stsp13_lat_vary, stsp13_lat_defined, + stsp13_fwhm_value, stsp13_fwhm_range, stsp13_fwhm_dstep, stsp13_fwhm_vary, stsp13_fwhm_defined, + stsp13_tcen_value, stsp13_tcen_range, stsp13_tcen_dstep, stsp13_tcen_vary, stsp13_tcen_defined, + stsp21_long_value, stsp21_long_range, stsp21_long_dstep, stsp21_long_vary, stsp21_long_defined, + stsp21_lat_value, stsp21_lat_range, stsp21_lat_dstep, stsp21_lat_vary, stsp21_lat_defined, + stsp21_fwhm_value, stsp21_fwhm_range, stsp21_fwhm_dstep, stsp21_fwhm_vary, stsp21_fwhm_defined, + stsp21_tcen_value, stsp21_tcen_range, stsp21_tcen_dstep, stsp21_tcen_vary, stsp21_tcen_defined, + stsp22_long_value, stsp22_long_range, stsp22_long_dstep, stsp22_long_vary, stsp22_long_defined, + stsp22_lat_value, stsp22_lat_range, stsp22_lat_dstep, stsp22_lat_vary, stsp22_lat_defined, + stsp22_fwhm_value, stsp22_fwhm_range, stsp22_fwhm_dstep, stsp22_fwhm_vary, stsp22_fwhm_defined, + stsp22_tcen_value, stsp22_tcen_range, stsp22_tcen_dstep, stsp22_tcen_vary, stsp22_tcen_defined, + uesp_long1_value, uesp_long1_range, uesp_long1_dstep, uesp_long1_vary, uesp_long1_defined, + uesp_long2_value, uesp_long2_range, uesp_long2_dstep, uesp_long2_vary, uesp_long2_defined, + uesp_lathw_value, uesp_lathw_range, uesp_lathw_dstep, uesp_lathw_vary, uesp_lathw_defined, + uesp_taper_value, uesp_taper_range, uesp_taper_dstep, uesp_taper_vary, uesp_taper_defined, + uesp_temp_value, uesp_temp_range, uesp_temp_dstep, uesp_temp_vary, uesp_temp_defined, + // disc + rdisc1_value, rdisc1_range, rdisc1_dstep, rdisc1_vary, rdisc1_defined, + rdisc2_value, rdisc2_range, rdisc2_dstep, rdisc2_vary, rdisc2_defined, + height_disc_value, height_disc_range, height_disc_dstep, height_disc_vary, height_disc_defined, + beta_disc_value, beta_disc_range, beta_disc_dstep, beta_disc_vary, beta_disc_defined, + temp_disc_value, temp_disc_range, temp_disc_dstep, temp_disc_vary, temp_disc_defined, + texp_disc_value, texp_disc_range, texp_disc_dstep, texp_disc_vary, texp_disc_defined, + lin_limb_disc_value, lin_limb_disc_range, lin_limb_disc_dstep, lin_limb_disc_vary, lin_limb_disc_defined, + quad_limb_disc_value, quad_limb_disc_range, quad_limb_disc_dstep, quad_limb_disc_vary, quad_limb_disc_defined, + temp_edge_value, temp_edge_range, temp_edge_dstep, temp_edge_vary, temp_edge_defined, + absorb_edge_value, absorb_edge_range, absorb_edge_dstep, absorb_edge_vary, absorb_edge_defined, + //Bright-spot + radius_spot_value, radius_spot_range, radius_spot_dstep, radius_spot_vary, radius_spot_defined, + length_spot_value, length_spot_range, length_spot_dstep, length_spot_vary, length_spot_defined, + height_spot_value, height_spot_range, height_spot_dstep, height_spot_vary, height_spot_defined, + expon_spot_value, expon_spot_range, expon_spot_dstep, expon_spot_vary, expon_spot_defined, + epow_spot_value, epow_spot_range, epow_spot_dstep, epow_spot_vary, epow_spot_defined, + angle_spot_value, angle_spot_range, angle_spot_dstep, angle_spot_vary, angle_spot_defined, + yaw_spot_value, yaw_spot_range, yaw_spot_dstep, yaw_spot_vary, yaw_spot_defined, + temp_spot_value, temp_spot_range, temp_spot_dstep, temp_spot_vary, temp_spot_defined, + tilt_spot_value, tilt_spot_range, tilt_spot_dstep, tilt_spot_vary, tilt_spot_defined, + cfrac_spot_value, cfrac_spot_range, cfrac_spot_dstep, cfrac_spot_vary, cfrac_spot_defined, + // Computational parameters + delta_phase, nlat1f, nlat2f, nlat1c, nlat2c, npole, + nlatfill, nlngfill, lfudge, llo, lhi, phase1, phase2, nrad, wavelength, + roche1, roche2, eclipse1, eclipse2, glens1, use_radii, + tperiod, gdark_bolom1, gdark_bolom2, mucrit1, mucrit2, + pslimb1, pslimb2, mirror, add_disc, opaque, add_spot, nspot, iscale + ); + + + + // Compute light curve + double wdwarf, chisq, wnok, logg1, logg2, rv1, rv2; + Lcurve::pylight_curve_comp(model, times, exposes, ndivs, Tsize, info, calc, + lcstar1, lcdisc, lcedge, lcspot, lcstar2, + wdwarf, logg1, logg2, rv1, rv2, parallel_threshold); + + wdwarflogrv[0] = wdwarf; + wdwarflogrv[1] = logg1; + wdwarflogrv[2] = logg2; + wdwarflogrv[3] = rv1; + wdwarflogrv[4] = rv2; + if(info){ + std::cout << "White dwarf's contribution = " << wdwarf + << std::endl; + std::cout << "log10(g1 [cgs]) = " << logg1 << std::endl; + std::cout << "log10(g2 [cgs]) = " << logg2 << std::endl; + std::cout << "Vol-averaged r1 = " << rv1 << std::endl; + std::cout << "Vol-averaged r2 = " << rv2 << std::endl; + } + + } //try + catch(const Roche::Roche_Error& err){ + std::cerr << "Roche::Roche_Error exception thrown" << std::endl; + std::cerr << "libpylcurve: " << err.what() << std::endl; + exit(EXIT_FAILURE); + } + catch(const Lcurve::Lcurve_Error& err){ + std::cerr << "Lcurve::Lcurve_Error exception thrown" << std::endl; + std::cerr << err << std::endl; + exit(EXIT_FAILURE); + } + catch(const std::string& err){ + std::cerr << err << std::endl; + exit(EXIT_FAILURE); + } + catch(...){ + std::cerr << "Unknown exception caught in lroche" << std::endl; + exit(EXIT_FAILURE); + } + }//void pylcurve( + + + void pylcurve_smodel(const char* psmodel, + double *times, double *exposes, int *ndivs, int Tsize, + bool info, + double *calc, double *lcstar1, double *lcdisc, + double *lcedge, double *lcspot, double *lcstar2, + double *wdwarflogrv, int parallel_threshold){ + try{ + const std::string smodel = psmodel; + Lcurve::Model model(smodel); + //free(psmodel); + // Compute light curve + double wdwarf, logg1, logg2, rv1, rv2; + Lcurve::pylight_curve_comp(model, times, exposes, ndivs, Tsize, info, calc, + lcstar1, lcdisc, lcedge, lcspot, lcstar2, + wdwarf, logg1, logg2, rv1, rv2, parallel_threshold); + wdwarflogrv[0] = wdwarf; + wdwarflogrv[1] = logg1; + wdwarflogrv[2] = logg2; + wdwarflogrv[3] = rv1; + wdwarflogrv[4] = rv2; + if(info){ + std::cout << "White dwarf's contribution = " << wdwarf + << std::endl; + std::cout << "log10(g1 [cgs]) = " << logg1 << std::endl; + std::cout << "log10(g2 [cgs]) = " << logg2 << std::endl; + std::cout << "Vol-averaged r1 = " << rv1 << std::endl; + std::cout << "Vol-averaged r2 = " << rv2 << std::endl; + } + + } //try + catch(const Roche::Roche_Error& err){ + std::cerr << "Roche::Roche_Error exception thrown" << std::endl; + std::cerr << "libpylcurve: " << err.what() << std::endl; + exit(EXIT_FAILURE); + } + catch(const Lcurve::Lcurve_Error& err){ + std::cerr << "Lcurve::Lcurve_Error exception thrown" << std::endl; + std::cerr << err << std::endl; + exit(EXIT_FAILURE); + } + catch(const std::string& err){ + std::cerr << err << std::endl; + exit(EXIT_FAILURE); + } + catch(...){ + std::cerr << "Unknown exception caught in lroche" << std::endl; + exit(EXIT_FAILURE); + } + }//pylcurve_file +} //extern "C" diff --git a/src/pylight_curve_comp.cc b/src/pylight_curve_comp.cc new file mode 100644 index 0000000..bb9b2c9 --- /dev/null +++ b/src/pylight_curve_comp.cc @@ -0,0 +1,312 @@ +#include "trm/lcurve.h" +#include "trm/buffer2d.h" +#include "trm/constants.h" + +#ifdef _OPENMP + #include +#endif + +/** This routine computes the light curve corresponding to a particular + * model and times defined by some data. Data with negative or zero errors + * are skipped for speed, corresponding fit values are set = 0. + * \param mdl the model + * \param data the data defining the times + * \param scale work out linear scale factors by svd or not. This will either be a single number for all components + * or a different one for each, depending upon the modle parameter iscale. + * \param info if true, it prints out array sizes to stderr + * \param calc the computed light curve + * \lcstar1 the light curve of star1 + * \lcdisc the light curve of star1 + * \lcedge the light curve of star1 + * \lcspot the light curve of star1 + * \lcstar2 the light curve of star1 + * \param wdwarf contribution of the white at phase 0.5 + * \param logg2 flux-weighted logg for star 2, CGS units + * \param rv1 and rv2 are volume-averaged radii + * This file is modified by lijiao + */ + +void Lcurve::pylight_curve_comp(const Lcurve::Model& mdl, + double *time, double *exposes, int *ndivs, int Tsize, + bool info, + double *calc, double *lcstar1, double *lcdisc, + double *lcedge, double *lcspot, double *lcstar2, + double& wdwarf, + double& logg1, double& logg2, double& rv1, double& rv2, int parallel_threshold){ + + double r1, r2; + mdl.get_r1r2(r1, r2); + double rl2 = 1.- Roche::xl12(mdl.q, mdl.spin2); + if(r2 < 0) + r2 = rl2; + else if(r2 > rl2) + throw Lcurve_Error("pylight_curve_comp: the secondary star is larger than its Roche lobe!"); + + LDC ldc1 = mdl.get_ldc1(); + LDC ldc2 = mdl.get_ldc2(); + + // Compute gravitational radius of star 1 if need be. An extra factor + // of 4 saves multiplication in the innermost loops later + double rlens1 = 0.; + if(mdl.glens1){ + // Compute G(M1+M2), SI, and the separation a, SI. + double gm = std::pow(1000.*mdl.velocity_scale,3)*mdl.tperiod*Constants::DAY/Constants::TWOPI; + double a = std::pow(gm/Subs::sqr(Constants::TWOPI/Constants::DAY/mdl.tperiod),1./3.); + rlens1 = 4.*gm/(1.+mdl.q)/a/Subs::sqr(Constants::C); + } + // std::cout<<"lijiao pylight_curve_comp rlens1 = " << rlens1< fstar2; + wdwarf = comp_star1(mdl.iangle, ldc1, 0.5, 0., 1, mdl.q, mdl.beam_factor1, + mdl.velocity_scale, gint, star1f, star1c); + + // calculate flux-weighted loggs and volume-averaged radii + logg1 = comp_gravity1(mdl, star1f); + logg2 = comp_gravity2(mdl, star2f); + rv1 = comp_radius1(star1f); + rv2 = comp_radius2(star2f); + + return; +} From cd96dff17cd3ee6559bc4bc0513819d40c5aaed8 Mon Sep 17 00:00:00 2001 From: lidihei <912219465@qq.com> Date: Wed, 19 Oct 2022 14:47:40 +0800 Subject: [PATCH 03/12] add libpylcurve.la in Makefile.am for python wrapper --- src/Makefile.am | 9 +++++++++ 1 file changed, 9 insertions(+) diff --git a/src/Makefile.am b/src/Makefile.am index 7844d50..51c3838 100644 --- a/src/Makefile.am +++ b/src/Makefile.am @@ -38,6 +38,15 @@ simplex_LDFLAGS = -fopenmp visualise_SOURCES = visualise.cc numface.cc set_star_grid.cc set_disc_grid.cc trm_lcurve.cc light_curve_comp.cc re_scale.cc set_star_continuum.cc \ comp_light.cc set_disc_continuum.cc set_bright_spot_grid.cc star_eclipse.cc comp_gravity.cc comp_radius.cc + +## create a library for python wrapper +lib_LTLIBRARIES = libpylcurve.la + +libpylcurve_la_SOURCES = libpylcurve.cc numface.cc set_star_grid.cc set_star_continuum.cc comp_light.cc trm_lcurve.cc pylight_curve_comp.cc re_scale.cc \ +set_disc_grid.cc set_disc_continuum.cc set_bright_spot_grid.cc star_eclipse.cc comp_gravity.cc comp_radius.cc light_curve_comp.cc +libpylcurve_la_CPPFLAGS = -fopenmp -shared -fPIC -I../include +libpylcurve_la_LDFLAGS = -fopenmp -shared -fPIC + AM_CPPFLAGS = -I../include DATE = $(shell date) From b1f15838cc8440812970428ef15f18afad406a8c Mon Sep 17 00:00:00 2001 From: lidihei <912219465@qq.com> Date: Wed, 19 Oct 2022 15:00:10 +0800 Subject: [PATCH 04/12] add pylight_curve_comp & Model Pparam from arguments to lcurve.h --- include/trm/lcurve.h | 109 +++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 109 insertions(+) diff --git a/include/trm/lcurve.h b/include/trm/lcurve.h index 1b69feb..3c55bdf 100644 --- a/include/trm/lcurve.h +++ b/include/trm/lcurve.h @@ -169,6 +169,15 @@ namespace Lcurve { if(!istr) defined = true; } + //! Constructor from arguments added by lijiao + Pparam(double value, double range, double dstep, bool vary, bool defined){ + this->value = value; + this->range = range; + this->dstep = dstep; + this->vary = vary; + this->defined = defined; + } + //! Implicit conversion to a double by just returning the value operator double () const { return value; } @@ -250,6 +259,97 @@ namespace Lcurve { //! Constructor from a file Model(const std::string& file); + //! Constructor from arguments added by lijiao + Model(//Binary and stars + double q_value, double q_range, double q_dstep, bool q_vary, bool q_defined, + double iangle_value, double iangle_range, double iangle_dstep, bool iangle_vary, bool iangle_defined, + double r1_value, double r1_range, double r1_dstep, bool r1_vary, bool r1_defined, + double r2_value, double r2_range, double r2_dstep, bool r2_vary, bool r2_defined, + double cphi3_value, double cphi3_range, double cphi3_dstep, bool cphi3_vary, bool cphi3_defined, + double cphi4_value, double cphi4_range, double cphi4_dstep, bool cphi4_vary, bool cphi4_defined, + double spin1_value, double spin1_range, double spin1_dstep, bool spin1_vary, bool spin1_defined, + double spin2_value, double spin2_range, double spin2_dstep, bool spin2_vary, bool spin2_defined, + double t1_value, double t1_range, double t1_dstep, bool t1_vary, bool t1_defined, + double t2_value, double t2_range, double t2_dstep, bool t2_vary, bool t2_defined, + double ldc1_1_value, double ldc1_1_range, double ldc1_1_dstep, bool ldc1_1_vary, bool ldc1_1_defined, + double ldc1_2_value, double ldc1_2_range, double ldc1_2_dstep, bool ldc1_2_vary, bool ldc1_2_defined, + double ldc1_3_value, double ldc1_3_range, double ldc1_3_dstep, bool ldc1_3_vary, bool ldc1_3_defined, + double ldc1_4_value, double ldc1_4_range, double ldc1_4_dstep, bool ldc1_4_vary, bool ldc1_4_defined, + double ldc2_1_value, double ldc2_1_range, double ldc2_1_dstep, bool ldc2_1_vary, bool ldc2_1_defined, + double ldc2_2_value, double ldc2_2_range, double ldc2_2_dstep, bool ldc2_2_vary, bool ldc2_2_defined, + double ldc2_3_value, double ldc2_3_range, double ldc2_3_dstep, bool ldc2_3_vary, bool ldc2_3_defined, + double ldc2_4_value, double ldc2_4_range, double ldc2_4_dstep, bool ldc2_4_vary, bool ldc2_4_defined, + double velocity_scale_value, double velocity_scale_range, double velocity_scale_dstep, bool velocity_scale_vary, bool velocity_scale_defined, + double beam_factor1_value, double beam_factor1_range, double beam_factor1_dstep, bool beam_factor1_vary, bool beam_factor1_defined, + double beam_factor2_value, double beam_factor2_range, double beam_factor2_dstep, bool beam_factor2_vary, bool beam_factor2_defined, + //General + double t0_value, double t0_range, double t0_dstep, bool t0_vary, bool t0_defined, + double period_value, double period_range, double period_dstep, bool period_vary, bool period_defined, + double pdot_value, double pdot_range, double pdot_dstep, bool pdot_vary, bool pdot_defined, + double deltat_value, double deltat_range, double deltat_dstep, bool deltat_vary, bool deltat_defined, + double gravity_dark1_value, double gravity_dark1_range, double gravity_dark1_dstep, bool gravity_dark1_vary, bool gravity_dark1_defined, + double gravity_dark2_value, double gravity_dark2_range, double gravity_dark2_dstep, bool gravity_dark2_vary, bool gravity_dark2_defined, + double absorb_value, double absorb_range, double absorb_dstep, bool absorb_vary, bool absorb_defined, + double slope_value, double slope_range, double slope_dstep, bool slope_vary, bool slope_defined, + double quad_value, double quad_range, double quad_dstep, bool quad_vary, bool quad_defined, + double cube_value, double cube_range, double cube_dstep, bool cube_vary, bool cube_defined, + double third_value, double third_range, double third_dstep, bool third_vary, bool third_defined, + // Star spots + double stsp11_long_value, double stsp11_long_range, double stsp11_long_dstep, bool stsp11_long_vary, bool stsp11_long_defined, + double stsp11_lat_value, double stsp11_lat_range, double stsp11_lat_dstep, bool stsp11_lat_vary, bool stsp11_lat_defined, + double stsp11_fwhm_value, double stsp11_fwhm_range, double stsp11_fwhm_dstep, bool stsp11_fwhm_vary, bool stsp11_fwhm_defined, + double stsp11_tcen_value, double stsp11_tcen_range, double stsp11_tcen_dstep, bool stsp11_tcen_vary, bool stsp11_tcen_defined, + double stsp12_long_value, double stsp12_long_range, double stsp12_long_dstep, bool stsp12_long_vary, bool stsp12_long_defined, + double stsp12_lat_value, double stsp12_lat_range, double stsp12_lat_dstep, bool stsp12_lat_vary, bool stsp12_lat_defined, + double stsp12_fwhm_value, double stsp12_fwhm_range, double stsp12_fwhm_dstep, bool stsp12_fwhm_vary, bool stsp12_fwhm_defined, + double stsp12_tcen_value, double stsp12_tcen_range, double stsp12_tcen_dstep, bool stsp12_tcen_vary, bool stsp12_tcen_defined, + double stsp13_long_value, double stsp13_long_range, double stsp13_long_dstep, bool stsp13_long_vary, bool stsp13_long_defined, + double stsp13_lat_value, double stsp13_lat_range, double stsp13_lat_dstep, bool stsp13_lat_vary, bool stsp13_lat_defined, + double stsp13_fwhm_value, double stsp13_fwhm_range, double stsp13_fwhm_dstep, bool stsp13_fwhm_vary, bool stsp13_fwhm_defined, + double stsp13_tcen_value, double stsp13_tcen_range, double stsp13_tcen_dstep, bool stsp13_tcen_vary, bool stsp13_tcen_defined, + double stsp21_long_value, double stsp21_long_range, double stsp21_long_dstep, bool stsp21_long_vary, bool stsp21_long_defined, + double stsp21_lat_value, double stsp21_lat_range, double stsp21_lat_dstep, bool stsp21_lat_vary, bool stsp21_lat_defined, + double stsp21_fwhm_value, double stsp21_fwhm_range, double stsp21_fwhm_dstep, bool stsp21_fwhm_vary, bool stsp21_fwhm_defined, + double stsp21_tcen_value, double stsp21_tcen_range, double stsp21_tcen_dstep, bool stsp21_tcen_vary, bool stsp21_tcen_defined, + double stsp22_long_value, double stsp22_long_range, double stsp22_long_dstep, bool stsp22_long_vary, bool stsp22_long_defined, + double stsp22_lat_value, double stsp22_lat_range, double stsp22_lat_dstep, bool stsp22_lat_vary, bool stsp22_lat_defined, + double stsp22_fwhm_value, double stsp22_fwhm_range, double stsp22_fwhm_dstep, bool stsp22_fwhm_vary, bool stsp22_fwhm_defined, + double stsp22_tcen_value, double stsp22_tcen_range, double stsp22_tcen_dstep, bool stsp22_tcen_vary, bool stsp22_tcen_defined, + double uesp_long1_value, double uesp_long1_range, double uesp_long1_dstep, bool uesp_long1_vary, bool uesp_long1_defined, + double uesp_long2_value, double uesp_long2_range, double uesp_long2_dstep, bool uesp_long2_vary, bool uesp_long2_defined, + double uesp_lathw_value, double uesp_lathw_range, double uesp_lathw_dstep, bool uesp_lathw_vary, bool uesp_lathw_defined, + double uesp_taper_value, double uesp_taper_range, double uesp_taper_dstep, bool uesp_taper_vary, bool uesp_taper_defined, + double uesp_temp_value, double uesp_temp_range, double uesp_temp_dstep, bool uesp_temp_vary, bool uesp_temp_defined, + // disc + double rdisc1_value, double rdisc1_range, double rdisc1_dstep, bool rdisc1_vary, bool rdisc1_defined, //disc + double rdisc2_value, double rdisc2_range, double rdisc2_dstep, bool rdisc2_vary, bool rdisc2_defined, + double height_disc_value, double height_disc_range, double height_disc_dstep, bool height_disc_vary, bool height_disc_defined, + double beta_disc_value, double beta_disc_range, double beta_disc_dstep, bool beta_disc_vary, bool beta_disc_defined, + double temp_disc_value, double temp_disc_range, double temp_disc_dstep, bool temp_disc_vary, bool temp_disc_defined, + double texp_disc_value, double texp_disc_range, double texp_disc_dstep, bool texp_disc_vary, bool texp_disc_defined, + double lin_limb_disc_value, double lin_limb_disc_range, double lin_limb_disc_dstep, bool lin_limb_disc_vary, bool lin_limb_disc_defined, + double quad_limb_disc_value, double quad_limb_disc_range, double quad_limb_disc_dstep, bool quad_limb_disc_vary, bool quad_limb_disc_defined, + double temp_edge_value, double temp_edge_range, double temp_edge_dstep, bool temp_edge_vary, bool temp_edge_defined, + double absorb_edge_value, double absorb_edge_range, double absorb_edge_dstep, bool absorb_edge_vary, bool absorb_edge_defined, + //Bright-spot + double radius_spot_value, double radius_spot_range, double radius_spot_dstep, bool radius_spot_vary, bool radius_spot_defined, //Bright-spot + double length_spot_value, double length_spot_range, double length_spot_dstep, bool length_spot_vary, bool length_spot_defined, + double height_spot_value, double height_spot_range, double height_spot_dstep, bool height_spot_vary, bool height_spot_defined, + double expon_spot_value, double expon_spot_range, double expon_spot_dstep, bool expon_spot_vary, bool expon_spot_defined, + double epow_spot_value, double epow_spot_range, double epow_spot_dstep, bool epow_spot_vary, bool epow_spot_defined, + double angle_spot_value, double angle_spot_range, double angle_spot_dstep, bool angle_spot_vary, bool angle_spot_defined, + double yaw_spot_value, double yaw_spot_range, double yaw_spot_dstep, bool yaw_spot_vary, bool yaw_spot_defined, + double temp_spot_value, double temp_spot_range, double temp_spot_dstep, bool temp_spot_vary, bool temp_spot_defined, + double tilt_spot_value, double tilt_spot_range, double tilt_spot_dstep, bool tilt_spot_vary, bool tilt_spot_defined, + double cfrac_spot_value, double cfrac_spot_range, double cfrac_spot_dstep, bool cfrac_spot_vary, bool cfrac_spot_defined, + // Computational parameters + double delta_phase, int nlat1f, int nlat2f, int nlat1c, int nlat2c, bool npole, + int nlatfill, int nlngfill, double lfudge, double llo, double lhi, double phase1, double phase2, int nrad, double wavelength, + bool roche1, bool roche2, bool eclipse1, bool eclipse2, bool glens1, bool use_radii, + double tperiod, double gdark_bolom1, double gdark_bolom2, double mucrit1, double mucrit2, + const char* pslimb1, const char* pslimb2, bool mirror, bool add_disc, bool opaque, bool add_spot, int nspot, bool iscale + ); + //! Number of variable parameters int nvary() const; @@ -857,6 +957,15 @@ namespace Lcurve { double& chisq, double& wnok, double& logg1, double& logg2, double& rv1, double& rv2); + //! Computes an entire light curve corresponding to a given time for python added by lijiao. + void pylight_curve_comp(const Lcurve::Model& mdl, + double *time, double *expose, int *ndiv, int Tsize, + bool info, + double *calc, double *lcstar1, double *lcdisc, + double *lcedge, double *lcspot, double *lcstar2, + double& wdwarf, + double& logg1, double& logg2, double& rv1, double& rv2, int parallel_threshold); + //! Re-scales a fit to minimise chi**2 double re_scale(const Lcurve::Data& data, Subs::Array1D& fit, double& chisq, double& wnok); From 744aed66a87abfdf26794c7bb41d39c155bb095f Mon Sep 17 00:00:00 2001 From: lidihei <912219465@qq.com> Date: Thu, 20 Oct 2022 21:48:13 +0800 Subject: [PATCH 05/12] add pyinstall --- pyinstall | 8 ++++++++ 1 file changed, 8 insertions(+) create mode 100755 pyinstall diff --git a/pyinstall b/pyinstall new file mode 100755 index 0000000..0979ec1 --- /dev/null +++ b/pyinstall @@ -0,0 +1,8 @@ +#!/bin/sh -v +# +#install the code for python wrapper of pylcurve (https://github.com/lidihei/pylcurve) + +autoreconf -f -i -s + +./configure --prefix=$TRM_SOFTWARE +make install From 58249101fa9848015702a6cde1b4fa8abc81a967 Mon Sep 17 00:00:00 2001 From: lidihei <912219465@qq.com> Date: Thu, 20 Oct 2022 22:32:53 +0800 Subject: [PATCH 06/12] libpylcurve.cc stsp21* --> stsp22 --- src/libpylcurve.cc | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/libpylcurve.cc b/src/libpylcurve.cc index c825d05..9272883 100644 --- a/src/libpylcurve.cc +++ b/src/libpylcurve.cc @@ -379,10 +379,10 @@ extern "C"{ double stsp21_lat_value, double stsp21_lat_range, double stsp21_lat_dstep, bool stsp21_lat_vary, bool stsp21_lat_defined, double stsp21_fwhm_value, double stsp21_fwhm_range, double stsp21_fwhm_dstep, bool stsp21_fwhm_vary, bool stsp21_fwhm_defined, double stsp21_tcen_value, double stsp21_tcen_range, double stsp21_tcen_dstep, bool stsp21_tcen_vary, bool stsp21_tcen_defined, - double stsp21_long_value, double stsp21_long_range, double stsp21_long_dstep, bool stsp21_long_vary, bool stsp21_long_defined, - double stsp21_lat_value, double stsp21_lat_range, double stsp21_lat_dstep, bool stsp21_lat_vary, bool stsp21_lat_defined, - double stsp21_fwhm_value, double stsp21_fwhm_range, double stsp21_fwhm_dstep, bool stsp21_fwhm_vary, bool stsp21_fwhm_defined, - double stsp21_tcen_value, double stsp21_tcen_range, double stsp21_tcen_dstep, bool stsp21_tcen_vary, bool stsp21_tcen_defined, + double stsp22_long_value, double stsp22_long_range, double stsp22_long_dstep, bool stsp22_long_vary, bool stsp22_long_defined, + double stsp22_lat_value, double stsp22_lat_range, double stsp22_lat_dstep, bool stsp22_lat_vary, bool stsp22_lat_defined, + double stsp22_fwhm_value, double stsp22_fwhm_range, double stsp22_fwhm_dstep, bool stsp22_fwhm_vary, bool stsp22_fwhm_defined, + double stsp22_tcen_value, double stsp22_tcen_range, double stsp22_tcen_dstep, bool stsp22_tcen_vary, bool stsp22_tcen_defined, double uesp_long1_value, double uesp_long1_range, double uesp_long1_dstep, bool uesp_long1_vary, bool uesp_long1_defined, double uesp_long2_value, double uesp_long2_range, double uesp_long2_dstep, bool uesp_long2_vary, bool uesp_long2_defined, double uesp_lathw_value, double uesp_lathw_range, double uesp_lathw_dstep, bool uesp_lathw_vary, bool uesp_lathw_defined, From 86c3c643116b9a49a615cc032aed08bc52d19e9a Mon Sep 17 00:00:00 2001 From: lidihei <912219465@qq.com> Date: Fri, 21 Oct 2022 11:07:31 +0800 Subject: [PATCH 07/12] python wrapper --- README.md | 18 ++++++++++++++++++ 1 file changed, 18 insertions(+) create mode 100644 README.md diff --git a/README.md b/README.md new file mode 100644 index 0000000..bd69db4 --- /dev/null +++ b/README.md @@ -0,0 +1,18 @@ + +'lcurve' is for modelling of light curves. + +Installation order: subs --> colly, binary --> roche --> lcurve + +The file called 'Lcurve' that is generated in this directory defines +aliases (csh-style) and prints a help message when sourced. + +Installation can be painful I'm afraid, but I don't have the time to make +a better job of it I am afraid. For details, please see: + + https://cygnus.astro.warwick.ac.uk/phsaap/software + + +This version allows python to call function in DLLs or shared libraries ([ctypes.CDLL](https://docs.python.org/3/library/ctypes.html)). +The python wrapper of lcurve is [pylcurve](https://github.com/lidihei/pylcurve), which is written by Jiao Li (lijiao@bao.ac.cn). + + From 7166a69717896a7cc2a3be8867d19ac8eaeba299 Mon Sep 17 00:00:00 2001 From: lidihei <912219465@qq.com> Date: Mon, 6 Feb 2023 17:01:28 +0800 Subject: [PATCH 08/12] gdark_bolom to bool --- README.md | 1 - include/trm/lcurve.h | 2 +- src/libpylcurve.cc | 2 +- src/trm_lcurve.cc | 2 +- 4 files changed, 3 insertions(+), 4 deletions(-) diff --git a/README.md b/README.md index bd69db4..0f92ae1 100644 --- a/README.md +++ b/README.md @@ -15,4 +15,3 @@ a better job of it I am afraid. For details, please see: This version allows python to call function in DLLs or shared libraries ([ctypes.CDLL](https://docs.python.org/3/library/ctypes.html)). The python wrapper of lcurve is [pylcurve](https://github.com/lidihei/pylcurve), which is written by Jiao Li (lijiao@bao.ac.cn). - diff --git a/include/trm/lcurve.h b/include/trm/lcurve.h index 3c55bdf..980e2ce 100644 --- a/include/trm/lcurve.h +++ b/include/trm/lcurve.h @@ -346,7 +346,7 @@ namespace Lcurve { double delta_phase, int nlat1f, int nlat2f, int nlat1c, int nlat2c, bool npole, int nlatfill, int nlngfill, double lfudge, double llo, double lhi, double phase1, double phase2, int nrad, double wavelength, bool roche1, bool roche2, bool eclipse1, bool eclipse2, bool glens1, bool use_radii, - double tperiod, double gdark_bolom1, double gdark_bolom2, double mucrit1, double mucrit2, + double tperiod, bool gdark_bolom1, bool gdark_bolom2, double mucrit1, double mucrit2, const char* pslimb1, const char* pslimb2, bool mirror, bool add_disc, bool opaque, bool add_spot, int nspot, bool iscale ); diff --git a/src/libpylcurve.cc b/src/libpylcurve.cc index 9272883..6171a2d 100644 --- a/src/libpylcurve.cc +++ b/src/libpylcurve.cc @@ -414,7 +414,7 @@ extern "C"{ double delta_phase, int nlat1f, int nlat2f, int nlat1c, int nlat2c, bool npole, int nlatfill, int nlngfill, double lfudge, double llo, double lhi, double phase1, double phase2, int nrad, double wavelength, bool roche1, bool roche2, bool eclipse1, bool eclipse2, bool glens1, bool use_radii, - double tperiod, double gdark_bolom1, double gdark_bolom2, double mucrit1, double mucrit2, + double tperiod, bool gdark_bolom1, bool gdark_bolom2, double mucrit1, double mucrit2, const char* pslimb1, const char* pslimb2, bool mirror, bool add_disc, bool opaque, bool add_spot, int nspot, bool iscale, bool info, int parallel_threshold ){ diff --git a/src/trm_lcurve.cc b/src/trm_lcurve.cc index 3cf5bfd..b9481a0 100644 --- a/src/trm_lcurve.cc +++ b/src/trm_lcurve.cc @@ -508,7 +508,7 @@ Lcurve::Model::Model(//Binary and stars double delta_phase, int nlat1f, int nlat2f, int nlat1c, int nlat2c, bool npole, int nlatfill, int nlngfill, double lfudge, double llo, double lhi, double phase1, double phase2, int nrad, double wavelength, bool roche1, bool roche2, bool eclipse1, bool eclipse2, bool glens1, bool use_radii, - double tperiod, double gdark_bolom1, double gdark_bolom2, double mucrit1, double mucrit2, + double tperiod, bool gdark_bolom1, bool gdark_bolom2, double mucrit1, double mucrit2, const char* pslimb1, const char* pslimb2, bool mirror, bool add_disc, bool opaque, bool add_spot, int nspot, bool iscale ) { From 8043a598bbfd42e060bc4fceb489a25eee68baae Mon Sep 17 00:00:00 2001 From: lidihei <912219465@qq.com> Date: Wed, 15 Feb 2023 10:36:12 +0800 Subject: [PATCH 09/12] check gdark_bolom if True GDCBOL = gravity_dark_value --- README | 4 ++-- src/set_star_continuum.cc | 2 +- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/README b/README index 24b1079..0f92ae1 100644 --- a/README +++ b/README @@ -12,6 +12,6 @@ a better job of it I am afraid. For details, please see: https://cygnus.astro.warwick.ac.uk/phsaap/software - - +This version allows python to call function in DLLs or shared libraries ([ctypes.CDLL](https://docs.python.org/3/library/ctypes.html)). +The python wrapper of lcurve is [pylcurve](https://github.com/lidihei/pylcurve), which is written by Jiao Li (lijiao@bao.ac.cn). diff --git a/src/set_star_continuum.cc b/src/set_star_continuum.cc index 886b785..69c9b10 100644 --- a/src/set_star_continuum.cc +++ b/src/set_star_continuum.cc @@ -54,7 +54,7 @@ void Lcurve::set_star_continuum(const Model& mdl, // used to give the desired behaviour of flux with gravity. const double GDCBOL1 = mdl.gdark_bolom1 ? mdl.gravity_dark1 : mdl.gravity_dark1 / Subs::dlpdlt(mdl.wavelength, mdl.t1); - + //std::cout<< "lijiao check GDCBOL1 = " << GDCBOL1 << std::endl; int nelem1 = star1.size(); // compute direction of star spot 11, 12, 13, and the uespot From 20078474d12f0266420c326a28c124088416fc4f Mon Sep 17 00:00:00 2001 From: lidihei <912219465@qq.com> Date: Fri, 17 Feb 2023 19:58:28 +0800 Subject: [PATCH 10/12] version description --- README.md | 10 ++++++++++ 1 file changed, 10 insertions(+) diff --git a/README.md b/README.md index 0f92ae1..9192d9f 100644 --- a/README.md +++ b/README.md @@ -1,3 +1,13 @@ +This is a modified version for pyclurve + +- major modification: +- - trm/lcurve.h add Pparam(double value, double range, double dstep, bool vary, bool defined) and + Constructor from arguments added by lijiao + Model(//Binary and stars + double q_value, double q_range, double q_dstep, bool q_vary, bool q_defined,... +- - src/Makefile.am add lib_LTLIBRARIES = libpylcurve.la ... +- - src/trm_lcurve.cc add Lcurve::Model::Model(//Binary and stars... +- - add libpylcurve.cc and pylight_curve_comp.cc 'lcurve' is for modelling of light curves. From fa1f6019ae480da12e6cb8d3e9fd8636d68b8a7d Mon Sep 17 00:00:00 2001 From: lidihei <912219465@qq.com> Date: Fri, 17 Feb 2023 20:00:30 +0800 Subject: [PATCH 11/12] version description --- README.md | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/README.md b/README.md index 9192d9f..2ec7f07 100644 --- a/README.md +++ b/README.md @@ -1,6 +1,6 @@ -This is a modified version for pyclurve +This is a modified version of lcurve for the python wrapper [pyclurve](https://github.com/lidihei/pylcurve) -- major modification: +- major modifications: - - trm/lcurve.h add Pparam(double value, double range, double dstep, bool vary, bool defined) and Constructor from arguments added by lijiao Model(//Binary and stars From 662d8c2df29be19675ca9d90915a12ea95ebb5a7 Mon Sep 17 00:00:00 2001 From: LiJiao Date: Thu, 13 Jun 2024 23:27:42 +0800 Subject: [PATCH 12/12] update_for_macOS --- .gitignore | 2 +- src/Makefile_Darwin.am | 61 ++++++++++++++++++++++++++++++++++++++++++ src/Makefile_Linux.am | 58 +++++++++++++++++++++++++++++++++++++++ 3 files changed, 120 insertions(+), 1 deletion(-) create mode 100644 src/Makefile_Darwin.am create mode 100644 src/Makefile_Linux.am diff --git a/.gitignore b/.gitignore index 1ca8787..f350931 100644 --- a/.gitignore +++ b/.gitignore @@ -3,7 +3,7 @@ Makefile Makefile.in *.o *.Po -src/* +#src/* !*.cc scripts/0651/* m4/* diff --git a/src/Makefile_Darwin.am b/src/Makefile_Darwin.am new file mode 100644 index 0000000..dcf11c7 --- /dev/null +++ b/src/Makefile_Darwin.am @@ -0,0 +1,61 @@ +## Process this file with automake to generate Makefile.in +## +## This is the file that must be edited if you are changing anything in the source directory + +## Programs + +progdir = @bindir@/@PACKAGE@ + +prog_PROGRAMS = lprofile visualise picture simplex lroche levmarq rotprof lroches + +lprofile_SOURCES = lprofile.cc numface.cc set_star_grid.cc set_star_emission.cc light_curve_comp.cc comp_light.cc trm_lcurve.cc \ +set_star_continuum.cc re_scale.cc set_disc_grid.cc set_bright_spot_grid.cc set_disc_continuum.cc star_eclipse.cc comp_gravity.cc \ +comp_radius.cc + +levmarq_SOURCES = levmarq.cc trm_lcurve.cc set_star_grid.cc comp_light.cc set_star_continuum.cc numface.cc light_curve_comp.cc re_scale.cc \ +set_disc_grid.cc set_disc_continuum.cc set_bright_spot_grid.cc star_eclipse.cc comp_gravity.cc comp_radius.cc +#levmarq_CPPFLAGS = -fopenmp -I../include +#levmarq_LDFLAGS = -fopenmp +levmarq_CPPFLAGS = -I../include + +lroche_SOURCES = lroche.cc numface.cc set_star_grid.cc set_star_continuum.cc comp_light.cc trm_lcurve.cc light_curve_comp.cc re_scale.cc \ +set_disc_grid.cc set_disc_continuum.cc set_bright_spot_grid.cc star_eclipse.cc comp_gravity.cc comp_radius.cc +#lroche_CPPFLAGS = -fopenmp -I../include +#lroche_LDFLAGS = -fopenmp +lroche_CPPFLAGS = -I../include + +lroches_SOURCES = lroche.cc numface.cc set_star_grid.cc set_star_continuum.cc comp_light.cc trm_lcurve.cc light_curve_comp.cc re_scale.cc \ +set_disc_grid.cc set_disc_continuum.cc set_bright_spot_grid.cc star_eclipse.cc comp_gravity.cc comp_radius.cc + +picture_SOURCES = picture.cc pos_disc.cc disc_eclipse.cc + +rotprof_SOURCES = rotprof.cc trm_lcurve.cc set_star_grid.cc light_curve_comp.cc re_scale.cc numface.cc set_star_continuum.cc \ +comp_light.cc set_disc_grid.cc star_eclipse.cc set_disc_continuum.cc set_bright_spot_grid.cc comp_gravity.cc comp_radius.cc + +simplex_SOURCES = simplex.cc trm_lcurve.cc set_star_grid.cc comp_light.cc set_star_continuum.cc numface.cc light_curve_comp.cc re_scale.cc \ +set_disc_grid.cc set_disc_continuum.cc set_bright_spot_grid.cc star_eclipse.cc comp_gravity.cc comp_radius.cc +#simplex_CPPFLAGS = -fopenmp -I../include +#simplex_LDFLAGS = -fopenmp +simplex_CPPFLAGS = -I../include +simplex_LDFLAGS = + +visualise_SOURCES = visualise.cc numface.cc set_star_grid.cc set_disc_grid.cc trm_lcurve.cc light_curve_comp.cc re_scale.cc set_star_continuum.cc \ +comp_light.cc set_disc_continuum.cc set_bright_spot_grid.cc star_eclipse.cc comp_gravity.cc comp_radius.cc + + +## create a library for python wrapper +lib_LTLIBRARIES = libpylcurve.la + +libpylcurve_la_SOURCES = libpylcurve.cc numface.cc set_star_grid.cc set_star_continuum.cc comp_light.cc trm_lcurve.cc pylight_curve_comp.cc re_scale.cc \ +set_disc_grid.cc set_disc_continuum.cc set_bright_spot_grid.cc star_eclipse.cc comp_gravity.cc comp_radius.cc light_curve_comp.cc +libpylcurve_la_CPPFLAGS = -fopenmp -shared -fPIC -I../include +libpylcurve_la_LDFLAGS = -fopenmp -shared -fPIC +libpylcurve_la_CPPFLAGS = -shared -fPIC -I../include +libpylcurve_la_LDFLAGS = -shared -fPIC + +AM_CPPFLAGS = -I../include + +DATE = $(shell date) + +install-data-hook: + echo "This is $(PACKAGE)-$(VERSION), built on $(DATE)" > $(progdir)/VERSION diff --git a/src/Makefile_Linux.am b/src/Makefile_Linux.am new file mode 100644 index 0000000..51c3838 --- /dev/null +++ b/src/Makefile_Linux.am @@ -0,0 +1,58 @@ +## Process this file with automake to generate Makefile.in +## +## This is the file that must be edited if you are changing anything in the source directory + +## Programs + +progdir = @bindir@/@PACKAGE@ + +prog_PROGRAMS = lprofile visualise picture simplex lroche levmarq rotprof lroches + +lprofile_SOURCES = lprofile.cc numface.cc set_star_grid.cc set_star_emission.cc light_curve_comp.cc comp_light.cc trm_lcurve.cc \ +set_star_continuum.cc re_scale.cc set_disc_grid.cc set_bright_spot_grid.cc set_disc_continuum.cc star_eclipse.cc comp_gravity.cc \ +comp_radius.cc + +levmarq_SOURCES = levmarq.cc trm_lcurve.cc set_star_grid.cc comp_light.cc set_star_continuum.cc numface.cc light_curve_comp.cc re_scale.cc \ +set_disc_grid.cc set_disc_continuum.cc set_bright_spot_grid.cc star_eclipse.cc comp_gravity.cc comp_radius.cc +levmarq_CPPFLAGS = -fopenmp -I../include +levmarq_LDFLAGS = -fopenmp + +lroche_SOURCES = lroche.cc numface.cc set_star_grid.cc set_star_continuum.cc comp_light.cc trm_lcurve.cc light_curve_comp.cc re_scale.cc \ +set_disc_grid.cc set_disc_continuum.cc set_bright_spot_grid.cc star_eclipse.cc comp_gravity.cc comp_radius.cc +lroche_CPPFLAGS = -fopenmp -I../include +lroche_LDFLAGS = -fopenmp + +lroches_SOURCES = lroche.cc numface.cc set_star_grid.cc set_star_continuum.cc comp_light.cc trm_lcurve.cc light_curve_comp.cc re_scale.cc \ +set_disc_grid.cc set_disc_continuum.cc set_bright_spot_grid.cc star_eclipse.cc comp_gravity.cc comp_radius.cc + +picture_SOURCES = picture.cc pos_disc.cc disc_eclipse.cc + +rotprof_SOURCES = rotprof.cc trm_lcurve.cc set_star_grid.cc light_curve_comp.cc re_scale.cc numface.cc set_star_continuum.cc \ +comp_light.cc set_disc_grid.cc star_eclipse.cc set_disc_continuum.cc set_bright_spot_grid.cc comp_gravity.cc comp_radius.cc + +simplex_SOURCES = simplex.cc trm_lcurve.cc set_star_grid.cc comp_light.cc set_star_continuum.cc numface.cc light_curve_comp.cc re_scale.cc \ +set_disc_grid.cc set_disc_continuum.cc set_bright_spot_grid.cc star_eclipse.cc comp_gravity.cc comp_radius.cc +simplex_CPPFLAGS = -fopenmp -I../include +simplex_LDFLAGS = -fopenmp + +visualise_SOURCES = visualise.cc numface.cc set_star_grid.cc set_disc_grid.cc trm_lcurve.cc light_curve_comp.cc re_scale.cc set_star_continuum.cc \ +comp_light.cc set_disc_continuum.cc set_bright_spot_grid.cc star_eclipse.cc comp_gravity.cc comp_radius.cc + + +## create a library for python wrapper +lib_LTLIBRARIES = libpylcurve.la + +libpylcurve_la_SOURCES = libpylcurve.cc numface.cc set_star_grid.cc set_star_continuum.cc comp_light.cc trm_lcurve.cc pylight_curve_comp.cc re_scale.cc \ +set_disc_grid.cc set_disc_continuum.cc set_bright_spot_grid.cc star_eclipse.cc comp_gravity.cc comp_radius.cc light_curve_comp.cc +libpylcurve_la_CPPFLAGS = -fopenmp -shared -fPIC -I../include +libpylcurve_la_LDFLAGS = -fopenmp -shared -fPIC + +AM_CPPFLAGS = -I../include + +DATE = $(shell date) + +install-data-hook: + echo "This is $(PACKAGE)-$(VERSION), built on $(DATE)" > $(progdir)/VERSION + + +