diff --git a/CMakeLists.txt b/CMakeLists.txt index 4761450..fa049c3 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -1,6 +1,6 @@ cmake_minimum_required(VERSION 3.15.2) # 设置项目名称与语言 -project(GCTL_POTENTIAL VERSION 1.0) +project(GCTL_POTENTIAL VERSION 2.0) # 添加配置配件编写的函数 include(CMakePackageConfigHelpers) diff --git a/lib/potential/gkernel_block.h b/lib/potential/gkernel_block.h index fa25a48..0c4f854 100644 --- a/lib/potential/gkernel_block.h +++ b/lib/potential/gkernel_block.h @@ -29,6 +29,9 @@ #define _GCTL_GRAV_KERNEL_BLOCK_H #include "gm_data.h" +#include "gctl/poly/block.h" +#include "gctl/utility/progress_bar.h" +#include "gctl/math/gmath.h" namespace gctl { diff --git a/lib/potential/gkernel_polygon.h b/lib/potential/gkernel_polygon.h index c58e32c..659c21a 100644 --- a/lib/potential/gkernel_polygon.h +++ b/lib/potential/gkernel_polygon.h @@ -29,6 +29,7 @@ #define _GCTL_GRAV_KERNEL_POLYGON_H #include "gm_data.h" +#include "gctl/poly/vertex.h" namespace gctl { diff --git a/lib/potential/gkernel_rectangle2d.h b/lib/potential/gkernel_rectangle2d.h index cb85157..9fec801 100644 --- a/lib/potential/gkernel_rectangle2d.h +++ b/lib/potential/gkernel_rectangle2d.h @@ -29,6 +29,9 @@ #define _GCTL_GRAV_KERNEL_RECTANGLE2D_H #include "gm_data.h" +#include "gctl/poly/point2c.h" +#include "gctl/poly/rectangle2d.h" +#include "gctl/utility/progress_bar.h" namespace gctl { diff --git a/lib/potential/gkernel_sphere.h b/lib/potential/gkernel_sphere.h index 256a99b..c19a6c7 100644 --- a/lib/potential/gkernel_sphere.h +++ b/lib/potential/gkernel_sphere.h @@ -29,6 +29,7 @@ #define _GCTL_GRAV_KERNEL_SPHERE_H #include "gm_data.h" +#include "gctl/poly/sphere.h" #include "gctl/optimization/loss_func.h" #ifdef GCTL_POTENTIAL_AUTODIFF diff --git a/lib/potential/gkernel_tesseroid.h b/lib/potential/gkernel_tesseroid.h index 459f151..b20081c 100644 --- a/lib/potential/gkernel_tesseroid.h +++ b/lib/potential/gkernel_tesseroid.h @@ -29,6 +29,8 @@ #define _GCTL_GRAV_KERNEL_TESSEROID_H #include "gm_data.h" +#include "gctl/poly/tesseroid.h" +#include "gctl/utility/progress_bar.h" #ifdef GCTL_POTENTIAL_TESS diff --git a/lib/potential/gkernel_tetrahedron.cpp b/lib/potential/gkernel_tetrahedron.cpp index 57edb22..2faf2dd 100644 --- a/lib/potential/gkernel_tetrahedron.cpp +++ b/lib/potential/gkernel_tetrahedron.cpp @@ -1039,7 +1039,7 @@ double gkernel_tetra_potential_sig_sph(const gctl::grav_tetrahedron &a_ele, cons gctl::gravtet_para* gp = a_ele.att; - op_c = a_op.s2c(); + op_c = s2c(a_op); int *v_order = a_ele.vec_order; @@ -1108,7 +1108,7 @@ gctl::point3dc gkernel_tetra_gradient_sig_sph(const gctl::grav_tetrahedron &a_el R[2][1] = cos((2.0+a_op.lon/180.0)*GCTL_Pi); R[2][2] = 0.0; - op_c = a_op.s2c(); + op_c = s2c(a_op); int *v_order = a_ele.vec_order; @@ -1178,7 +1178,7 @@ gctl::tensor gkernel_tetra_tensor_sig_sph(const gctl::grav_tetrahedron &a_ele, c R_T = R.transpose(); - op_c = a_op.s2c(); + op_c = s2c(a_op); int *v_order = a_ele.vec_order; diff --git a/lib/potential/gkernel_tetrahedron.h b/lib/potential/gkernel_tetrahedron.h index 65ad1de..dbdc64a 100644 --- a/lib/potential/gkernel_tetrahedron.h +++ b/lib/potential/gkernel_tetrahedron.h @@ -29,6 +29,9 @@ #define _GCTL_GRAV_KERNEL_TETRAHEDRON_H #include "gm_data.h" +#include "gctl/poly/tetrahedron.h" +#include "gctl/math/geometry3d.h" +#include "gctl/utility/progress_bar.h" namespace gctl { diff --git a/lib/potential/gkernel_triangle.cpp b/lib/potential/gkernel_triangle.cpp index d0e5034..361ac26 100644 --- a/lib/potential/gkernel_triangle.cpp +++ b/lib/potential/gkernel_triangle.cpp @@ -360,7 +360,7 @@ double gkernel_triangle_potential_sig_sph(const gctl::grav_triangle &a_ele, cons gctl::gravtri_para *gp = a_ele.att; - op_c = a_op.s2c(); + op_c = s2c(a_op); r_ijk[0] = *a_ele.vert[0] - op_c; r_ijk[1] = *a_ele.vert[1] - op_c; @@ -423,7 +423,7 @@ gctl::point3dc gkernel_triangle_gradient_sig_sph(const gctl::grav_triangle &a_el R[2][1] = cos((2.0+a_op.lon/180.0)*GCTL_Pi); R[2][2] = 0.0; - op_c = a_op.s2c(); + op_c = s2c(a_op); r_ijk[0] = *a_ele.vert[0] - op_c; r_ijk[1] = *a_ele.vert[1] - op_c; @@ -487,7 +487,7 @@ gctl::tensor gkernel_triangle_tensor_sig_sph(const gctl::grav_triangle &a_ele, c R_T = R.transpose(); - op_c = a_op.s2c(); + op_c = s2c(a_op); r_ijk[0] = *a_ele.vert[0] - op_c; r_ijk[1] = *a_ele.vert[1] - op_c; diff --git a/lib/potential/gkernel_triangle.h b/lib/potential/gkernel_triangle.h index 51e12ee..048b1cd 100644 --- a/lib/potential/gkernel_triangle.h +++ b/lib/potential/gkernel_triangle.h @@ -29,6 +29,8 @@ #define _GCTL_GRAV_KERNEL_TRIANGLE_H #include "gm_data.h" +#include "gctl/poly/triangle.h" +#include "gctl/utility/progress_bar.h" namespace gctl { diff --git a/lib/potential/gkernel_triangle2d.cpp b/lib/potential/gkernel_triangle2d.cpp index bda0f9d..af184a4 100644 --- a/lib/potential/gkernel_triangle2d.cpp +++ b/lib/potential/gkernel_triangle2d.cpp @@ -443,7 +443,7 @@ void gkernel_triangle2d_vz2(gctl::matrix &out_kernel, const gctl::array< if (verbose == gctl::FullMsg) bar.progressed(i); else if (verbose == gctl::ShortMsg) bar.progressed_simple(i); - obsc = obsp[i].p2c(); + obsc = p2c(obsp[i]); #pragma omp parallel for private (j, obsg) schedule(guided) for (j = 0; j < e_size; j++) @@ -471,7 +471,7 @@ void gkernel_triangle2d_vx2(gctl::matrix &out_kernel, const gctl::array< if (verbose == gctl::FullMsg) bar.progressed(i); else if (verbose == gctl::ShortMsg) bar.progressed_simple(i); - obsc = obsp[i].p2c(); + obsc = p2c(obsp[i]); #pragma omp parallel for private (j, obsg) schedule(guided) for (j = 0; j < e_size; j++) @@ -500,7 +500,7 @@ void gkernel_triangle2d_vxx2(gctl::matrix &out_kernel, const gctl::array if (verbose == gctl::FullMsg) bar.progressed(i); else if (verbose == gctl::ShortMsg) bar.progressed_simple(i); - obsc = obsp[i].p2c(); + obsc = p2c(obsp[i]); #pragma omp parallel for private (j, t) schedule(guided) for (j = 0; j < e_size; j++) @@ -532,7 +532,7 @@ void gkernel_triangle2d_vxz2(gctl::matrix &out_kernel, const gctl::array if (verbose == gctl::FullMsg) bar.progressed(i); else if (verbose == gctl::ShortMsg) bar.progressed_simple(i); - obsc = obsp[i].p2c(); + obsc = p2c(obsp[i]); #pragma omp parallel for private (j, t) schedule(guided) for (j = 0; j < e_size; j++) @@ -564,7 +564,7 @@ void gkernel_triangle2d_vzx2(gctl::matrix &out_kernel, const gctl::array if (verbose == gctl::FullMsg) bar.progressed(i); else if (verbose == gctl::ShortMsg) bar.progressed_simple(i); - obsc = obsp[i].p2c(); + obsc = p2c(obsp[i]); #pragma omp parallel for private (j, t) schedule(guided) for (j = 0; j < e_size; j++) @@ -596,7 +596,7 @@ void gkernel_triangle2d_vzz2(gctl::matrix &out_kernel, const gctl::array if (verbose == gctl::FullMsg) bar.progressed(i); else if (verbose == gctl::ShortMsg) bar.progressed_simple(i); - obsc = obsp[i].p2c(); + obsc = p2c(obsp[i]); #pragma omp parallel for private (j, t) schedule(guided) for (j = 0; j < e_size; j++) @@ -774,7 +774,7 @@ void gobser_triangle2d_vz2(gctl::array &out_obs, const gctl::array &out_obs, const gctl::array &out_obs, const gctl::array &out_obs, const gctl::array &out_obs, const gctl::array &out_obs, const gctl::array &out_mGrad, #pragma omp parallel for private (j) shared(i) schedule(guided) for (j = 0; j < v_size; j++) { - if (angle(ops[i].s2c(), verts_ptr->at(j)) <= ang_limit*GCTL_Pi/180.0) + if (angle(s2c(ops[i]), verts_ptr->at(j)) <= ang_limit*GCTL_Pi/180.0) out_mGrad.at(j) += obs_diff[i] * gobser_mGrad_tricone_vr_sig(vert_neighList[j], j, ops[i], rho); } } @@ -379,7 +379,7 @@ void gobser_mGrad_tricone_vrp(gctl::array &out_mGrad, #pragma omp parallel for private (j) shared(i) schedule(guided) for (j = 0; j < v_size; j++) { - if (angle(ops[i].s2c(), verts_ptr->at(j)) <= ang_limit*GCTL_Pi/180.0) + if (angle(s2c(ops[i]), verts_ptr->at(j)) <= ang_limit*GCTL_Pi/180.0) out_mGrad.at(j) += obs_diff[i] * gobser_mGrad_tricone_vrp_sig(vert_neighList[j], j, ops[i], rho); } } @@ -424,7 +424,7 @@ void gobser_mGrad_tricone_vrt(gctl::array &out_mGrad, #pragma omp parallel for private (j) shared(i) schedule(guided) for (j = 0; j < v_size; j++) { - if (angle(ops[i].s2c(), verts_ptr->at(j)) <= ang_limit*GCTL_Pi/180.0) + if (angle(s2c(ops[i]), verts_ptr->at(j)) <= ang_limit*GCTL_Pi/180.0) out_mGrad.at(j) += obs_diff[i] * gobser_mGrad_tricone_vrt_sig(vert_neighList[j], j, ops[i], rho); } } @@ -469,7 +469,7 @@ void gobser_mGrad_tricone_vrr(gctl::array &out_mGrad, #pragma omp parallel for private (j) shared(i) schedule(guided) for (j = 0; j < v_size; j++) { - if (angle(ops[i].s2c(), verts_ptr->at(j)) <= ang_limit*GCTL_Pi/180.0) + if (angle(s2c(ops[i]), verts_ptr->at(j)) <= ang_limit*GCTL_Pi/180.0) out_mGrad.at(j) += obs_diff[i] * gobser_mGrad_tricone_vrr_sig(vert_neighList[j], j, ops[i], rho); } } @@ -631,7 +631,7 @@ double gobser_mGrad_tricone_vr_sig(const std::vector &a_lis R.y = sin((0.5-a_op.lat/180.0)*GCTL_Pi)*sin((2.0+a_op.lon/180.0)*GCTL_Pi); R.z = cos((0.5-a_op.lat/180.0)*GCTL_Pi); - op_c = a_op.s2c(); + op_c = s2c(a_op); double point_sum = 0.0; for (int j = 0; j < a_list.size(); j++) @@ -842,7 +842,7 @@ double gobser_mGrad_tricone_vrp_sig(const std::vector &a_li R.y = cos((2.0+a_op.lon/180.0)*GCTL_Pi); R.z = 0.0; - op_c = a_op.s2c(); + op_c = s2c(a_op); double point_sum = 0.0; for (int j = 0; j < a_list.size(); j++) @@ -1033,7 +1033,7 @@ double gobser_mGrad_tricone_vrt_sig(const std::vector &a_li R.y = cos((0.5-a_op.lat/180.0)*GCTL_Pi)*sin((2.0+a_op.lon/180.0)*GCTL_Pi); R.z = -1.0*sin((0.5-a_op.lat/180.0)*GCTL_Pi); - op_c = a_op.s2c(); + op_c = s2c(a_op); double point_sum = 0.0; for (int j = 0; j < a_list.size(); j++) @@ -1220,7 +1220,7 @@ double gobser_mGrad_tricone_vrr_sig(const std::vector &a_li R.y = sin((0.5-a_op.lat/180.0)*GCTL_Pi)*sin((2.0+a_op.lon/180.0)*GCTL_Pi); R.z = cos((0.5-a_op.lat/180.0)*GCTL_Pi); - op_c = a_op.s2c(); + op_c = s2c(a_op); double point_sum = 0.0; for (int j = 0; j < a_list.size(); j++) @@ -1404,7 +1404,7 @@ double gkernel_tricone_gji_pot_sig(const gctl::gravcone_gji &a_ele, const gctl:: gctl::gravcone_para_gji* gp = a_ele.att; - op_c = a_op.s2c(); + op_c = s2c(a_op); face_sum = edge_sum = 0.0; for (f = 0; f < 4; f++) @@ -1461,7 +1461,7 @@ double gkernel_tricone_gji_vr_sig(const gctl::gravcone_gji &a_ele, const gctl::p R.y = sin((0.5-a_op.lat/180.0)*GCTL_Pi)*sin((2.0+a_op.lon/180.0)*GCTL_Pi); R.z = cos((0.5-a_op.lat/180.0)*GCTL_Pi); - op_c = a_op.s2c(); + op_c = s2c(a_op); face_sum = edge_sum = 0.0; for (f = 0; f < 4; f++) @@ -1525,7 +1525,7 @@ double gkernel_tricone_gji_vrp_sig(const gctl::gravcone_gji &a_ele, const gctl:: R[2][1] = cos((2.0+a_op.lon/180.0)*GCTL_Pi); R[2][2] = 0.0; - op_c = a_op.s2c(); + op_c = s2c(a_op); face_sum = edge_sum = 0.0; for (f = 0; f < 4; f++) @@ -1594,7 +1594,7 @@ double gkernel_tricone_gji_vrt_sig(const gctl::gravcone_gji &a_ele, const gctl:: R[2][1] = cos((2.0+a_op.lon/180.0)*GCTL_Pi); R[2][2] = 0.0; - op_c = a_op.s2c(); + op_c = s2c(a_op); face_sum = edge_sum = 0.0; for (f = 0; f < 4; f++) @@ -1663,7 +1663,7 @@ double gkernel_tricone_gji_vrr_sig(const gctl::gravcone_gji &a_ele, const gctl:: R[2][1] = cos((2.0+a_op.lon/180.0)*GCTL_Pi); R[2][2] = 0.0; - op_c = a_op.s2c(); + op_c = s2c(a_op); face_sum = edge_sum = 0.0; for (f = 0; f < 4; f++) diff --git a/lib/potential/gm_data.h b/lib/potential/gm_data.h index 92fa6f1..0116011 100644 --- a/lib/potential/gm_data.h +++ b/lib/potential/gm_data.h @@ -34,12 +34,12 @@ #include #endif -#include "gctl/core.h" -#include "gctl/utility.h" -#include "gctl/geometry.h" -#include "gctl/maths.h" -#include "gctl/algorithms.h" -#include "gctl/io.h" +#include "gctl/core/array.h" +#include "gctl/poly/point3s.h" +#include "gctl/poly/tensor.h" +#include "gctl/math/gmath.h" +#include "gctl/io/dsv_io.h" +#include "gctl/utility/utc_time.h" namespace gctl { diff --git a/lib/potential/gm_regular_grid.h b/lib/potential/gm_regular_grid.h index 4287db6..46306db 100644 --- a/lib/potential/gm_regular_grid.h +++ b/lib/potential/gm_regular_grid.h @@ -30,6 +30,8 @@ #include "gctl/mesh/regular_grid.h" #include "gm_data.h" +#include "gctl/math/extrapolate.h" +#include "gctl/math/space_filter.h" namespace gctl { diff --git a/lib/potential/mkernel_dipole.cpp b/lib/potential/mkernel_dipole.cpp index b58a3cd..8070a78 100644 --- a/lib/potential/mkernel_dipole.cpp +++ b/lib/potential/mkernel_dipole.cpp @@ -48,7 +48,7 @@ gctl::point3dc gctl::magkernel_single(const mag_dipole &a_dipole, const point3ds else R = transform_matrix(a_op); double r = a_op.rad; - point3dc pc = a_op.s2c(); + point3dc pc = s2c(a_op); // T -> nT 1e+9 * 1e-7 -> 1e+2 point3dc b = 1e+2*a_dipole.M*(3.0*dot(a_dipole.n, pc)/power5(r)*pc - 1.0/power3(r)*a_dipole.n); diff --git a/lib/potential/mkernel_dipole.h b/lib/potential/mkernel_dipole.h index a73e7b4..90604c7 100644 --- a/lib/potential/mkernel_dipole.h +++ b/lib/potential/mkernel_dipole.h @@ -29,6 +29,8 @@ #define _GCTL_MAG_KERNEL_DIPOLE_H #include "gm_data.h" +#include "gctl/poly/vertex.h" +#include "gctl/utility/progress_bar.h" namespace gctl { diff --git a/lib/potential/mkernel_tesseroid.h b/lib/potential/mkernel_tesseroid.h index 4112ef5..5092550 100644 --- a/lib/potential/mkernel_tesseroid.h +++ b/lib/potential/mkernel_tesseroid.h @@ -29,6 +29,8 @@ #define _GCTL_MAG_KERNEL_TESSEROID_H #include "gm_data.h" +#include "gctl/poly/tesseroid.h" +#include "gctl/utility/progress_bar.h" #ifdef GCTL_POTENTIAL_MAGTESS diff --git a/lib/potential/mkernel_tetrahedron.cpp b/lib/potential/mkernel_tetrahedron.cpp index c61c79b..1d6d0d6 100644 --- a/lib/potential/mkernel_tetrahedron.cpp +++ b/lib/potential/mkernel_tetrahedron.cpp @@ -142,7 +142,7 @@ void gctl::callink_magnetic_para_earth_sph(array &in_tet, array mag_z.y = cos(I)*cos(A); mag_z.z = -1.0*sin(I); - s = in_tet[i].center().c2s(); + s = c2s(in_tet[i].center()); // magnetic susceptibility is taken as one here // rotate the local coordinate system to the regular status out_para[i].B = field_tense * mag_z.rotate((90.0 - s.lat)*GCTL_Pi/180.0, 0.0, (90.0 + s.lon)*GCTL_Pi/180.0).normal(); @@ -231,7 +231,7 @@ gctl::point3dc gctl::magkernel_single(const gctl::mag_tetrahedron &a_ele, const else R = transform_matrix(a_op); gctl::magtet_para* gp = a_ele.att; - gctl::point3dc pc = a_op.s2c(); + gctl::point3dc pc = s2c(a_op); int *v_order = a_ele.vec_order; for (f = 0; f < 4; f++) diff --git a/lib/potential/mkernel_tetrahedron.h b/lib/potential/mkernel_tetrahedron.h index 90c7e2b..4886fca 100644 --- a/lib/potential/mkernel_tetrahedron.h +++ b/lib/potential/mkernel_tetrahedron.h @@ -29,6 +29,8 @@ #define _GCTL_MAG_KERNEL_TETRAHEDRON_H #include "gm_data.h" +#include "gctl/poly/tetrahedron.h" +#include "gctl/utility/progress_bar.h" namespace gctl { diff --git a/lib/potential/mkernel_tetrahedron_Ren2017.cpp b/lib/potential/mkernel_tetrahedron_Ren2017.cpp index e9c9774..d8982a1 100644 --- a/lib/potential/mkernel_tetrahedron_Ren2017.cpp +++ b/lib/potential/mkernel_tetrahedron_Ren2017.cpp @@ -148,7 +148,7 @@ void gctl::callink_magnetic_para_earth_sph(array &in_tet, array &out_kernel, const array &el #pragma omp parallel for private (j, obsp_c) schedule(guided) for (j = 0; j < e_size; j++) { - obsp_c = obsp[i].s2c(); + obsp_c = s2c(obsp[i]); out_kernel[i][j] = magkernel_tetrahedron_gradient_sig(ele[j], obsp_c); } } @@ -411,7 +411,7 @@ void gctl::magkernel(matrix &out_kernel, const array &ele, #pragma omp parallel for private (j, obsp_c) schedule(guided) for (j = 0; j < e_size; j++) { - obsp_c = obsp[i].s2c(); + obsp_c = s2c(obsp[i]); out_kernel[i][j] = magkernel_tetrahedron_tensor_sig(ele[j], obsp_c); } } @@ -461,7 +461,7 @@ gctl::point3dc gctl::magkernel_single(const magtet_ren17 &ele, const point3ds &o R[2][1] = cos((2.0+obsp.lon/180.0)*M_PI); R[2][2] = 0.0; - point3dc obsp_c = obsp.s2c(); + point3dc obsp_c = s2c(obsp); point3dc out_k = magkernel_tetrahedron_gradient_sig(ele, obsp_c); out_k = R * out_k; @@ -634,7 +634,7 @@ void gctl::magobser(array &out_obs, const array &ele, co #pragma omp parallel for private (i, obsp_c) schedule(guided) for (i = 0; i < o_size; i++) { - obsp_c = obsp[i].s2c(); + obsp_c = s2c(obsp[i]); out_obs[i] = out_obs[i] + sus[j] * magkernel_tetrahedron_gradient_sig(ele[j], obsp_c); } } @@ -680,7 +680,7 @@ void gctl::magobser(array &out_obs, const array &ele, cons #pragma omp parallel for private (i, obsp_c) schedule(guided) for (i = 0; i < o_size; i++) { - obsp_c = obsp[i].s2c(); + obsp_c = s2c(obsp[i]); out_obs[i] = out_obs[i] + sus[j] * magkernel_tetrahedron_tensor_sig(ele[j], obsp_c); } } diff --git a/lib/potential/mkernel_tetrahedron_Ren2017.h b/lib/potential/mkernel_tetrahedron_Ren2017.h index e08307a..4b33b60 100644 --- a/lib/potential/mkernel_tetrahedron_Ren2017.h +++ b/lib/potential/mkernel_tetrahedron_Ren2017.h @@ -29,6 +29,8 @@ #define _GCTL_MAG_KERNEL_TETRAHEDRON_REN2017_H #include "gm_data.h" +#include "gctl/poly/tetrahedron.h" +#include "gctl/utility/progress_bar.h" namespace gctl { diff --git a/lib/potential/mkernel_triangle.cpp b/lib/potential/mkernel_triangle.cpp index f0cbd49..7af12cb 100644 --- a/lib/potential/mkernel_triangle.cpp +++ b/lib/potential/mkernel_triangle.cpp @@ -137,7 +137,7 @@ void gctl::callink_magnetic_para_earth_sph(array &in_tri, array &in_tet, array &kernel, const array &top_ele, const for (j = 0; j < e_size; j++) { cen = 1.0/3.0*(*top_ele[j].vert[0] + *top_ele[j].vert[1] + *top_ele[j].vert[2]); - if (geometry3d::angle(obsp[i].s2c(), cen) < cut_angle*GCTL_Pi/180.0) + if (geometry3d::angle(s2c(obsp[i]), cen) < cut_angle*GCTL_Pi/180.0) { tmp_plt.r_id = i; tmp_plt.c_id = j; @@ -300,7 +300,7 @@ void gctl::magkernel(spmat &kernel, const array &top_ele, const for (j = 0; j < e_size; j++) { cen = 1.0/3.0*(*top_ele[j].vert[0] + *top_ele[j].vert[1] + *top_ele[j].vert[2]); - if (geometry3d::angle(obsp[i].s2c(), cen) < cut_angle*GCTL_Pi/180.0) + if (geometry3d::angle(s2c(obsp[i]), cen) < cut_angle*GCTL_Pi/180.0) { mag_b = magkernel_single(top_ele[j], obsp[i], Rs.get(i)); if (comp_type == Bx) triplts[tri_idx[j + i*e_size]].val = mag_b.y; @@ -319,7 +319,7 @@ void gctl::magkernel(spmat &kernel, const array &top_ele, const for (j = 0; j < e_size; j++) { cen = 1.0/3.0*(*btm_ele[j].vert[0] + *btm_ele[j].vert[1] + *btm_ele[j].vert[2]); - if (geometry3d::angle(obsp[i].s2c(), cen) < cut_angle*GCTL_Pi/180.0) + if (geometry3d::angle(s2c(obsp[i]), cen) < cut_angle*GCTL_Pi/180.0) { mag_b = magkernel_single(btm_ele[j], obsp[i], Rs.get(i)); if (comp_type == Bx) triplts[tri_idx[j + i*e_size]].val -= mag_b.y; diff --git a/lib/potential/mkernel_tricone.h b/lib/potential/mkernel_tricone.h index 1d99651..33036ba 100644 --- a/lib/potential/mkernel_tricone.h +++ b/lib/potential/mkernel_tricone.h @@ -29,6 +29,9 @@ #define _GCTL_MAG_KERNEL_TRICONE_H #include "gm_data.h" +#include "gctl/poly/tri_cone.h" +#include "gctl/utility/progress_bar.h" +#include "gctl/math/geometry3d.h" namespace gctl { diff --git a/lib/potential/mkernel_tricone_Ren2017.cpp b/lib/potential/mkernel_tricone_Ren2017.cpp index bccc695..27ce86c 100644 --- a/lib/potential/mkernel_tricone_Ren2017.cpp +++ b/lib/potential/mkernel_tricone_Ren2017.cpp @@ -107,7 +107,7 @@ void gctl::callink_magnetic_para_earth_sph(array &in_cone, array< mag_z.z = -1.0*sin(I); c = 1.0/3.0*(*in_cone[i].vert[0] + *in_cone[i].vert[1] + *in_cone[i].vert[2]); - s = c.c2s(); + s = c2s(c); // magnetic susceptibility is taken as one here // rotate the local coordinate system to the regular status mag_z = field_tense * mag_z.rotate((90.0 - s.lat)*M_PI/180.0, 0.0, (90.0 + s.lon)*M_PI/180.0).normal(); @@ -151,7 +151,7 @@ gctl::point3dc gctl::magkernel_single(const magcone_ren17 &a_ele, const point3ds else R = transform_matrix(a_op); // get the observation site in the local coordinates - gctl::point3dc site = a_op.s2c(); + gctl::point3dc site = s2c(a_op); gctl::point3dc out_grad(0.0, 0.0, 0.0); for (int f = 0; f < 4; ++f) { @@ -238,7 +238,7 @@ gctl::tensor gctl::magkernel_single_tensor(const magcone_ren17 &a_ele, const poi R_T = R.transpose(); // get the observation site in the local coordinates - gctl::point3dc site = a_op.s2c(); + gctl::point3dc site = s2c(a_op); gctl::tensor out_tensor(0.0); for (int f = 0; f < 4; ++f) { @@ -435,7 +435,7 @@ void gctl::magkernel(spmat &kernel, const array &top_ele, for (j = 0; j < e_size; j++) { cen = 1.0/3.0*(*top_ele[j].vert[0] + *top_ele[j].vert[1] + *top_ele[j].vert[2]); - if (geometry3d::angle(obsp[i].s2c(), cen) < cut_angle*GCTL_Pi/180.0) + if (geometry3d::angle(s2c(obsp[i]), cen) < cut_angle*GCTL_Pi/180.0) { tmp_plt.r_id = i; tmp_plt.c_id = j; @@ -457,7 +457,7 @@ void gctl::magkernel(spmat &kernel, const array &top_ele, for (j = 0; j < e_size; j++) { cen = 1.0/3.0*(*top_ele[j].vert[0] + *top_ele[j].vert[1] + *top_ele[j].vert[2]); - if (geometry3d::angle(obsp[i].s2c(), cen) < cut_angle*GCTL_Pi/180.0) + if (geometry3d::angle(s2c(obsp[i]), cen) < cut_angle*GCTL_Pi/180.0) { mag_b = magkernel_single(top_ele[j], obsp[i], Rs.get(i)); if (comp_type == Bx) triplts[tri_idx[j + i*e_size]].val = mag_b.y; @@ -476,7 +476,7 @@ void gctl::magkernel(spmat &kernel, const array &top_ele, for (j = 0; j < e_size; j++) { cen = 1.0/3.0*(*btm_ele[j].vert[0] + *btm_ele[j].vert[1] + *btm_ele[j].vert[2]); - if (geometry3d::angle(obsp[i].s2c(), cen) < cut_angle*GCTL_Pi/180.0) + if (geometry3d::angle(s2c(obsp[i]), cen) < cut_angle*GCTL_Pi/180.0) { mag_b = magkernel_single(btm_ele[j], obsp[i], Rs.get(i)); if (comp_type == Bx) triplts[tri_idx[j + i*e_size]].val -= mag_b.y; diff --git a/lib/potential/mkernel_tricone_Ren2017.h b/lib/potential/mkernel_tricone_Ren2017.h index e5595d7..7ea51d6 100644 --- a/lib/potential/mkernel_tricone_Ren2017.h +++ b/lib/potential/mkernel_tricone_Ren2017.h @@ -29,6 +29,9 @@ #define _GCTL_MAG_KERNEL_TRICONE_REN2017_H #include "gm_data.h" +#include "gctl/poly/tri_cone.h" +#include "gctl/math/geometry3d.h" +#include "gctl/utility/progress_bar.h" namespace gctl { diff --git a/tool/gridmanager-potential/gridmanager-potential.cpp b/tool/gridmanager-potential/gridmanager-potential.cpp index 2d9b9fe..246cb16 100644 --- a/tool/gridmanager-potential/gridmanager-potential.cpp +++ b/tool/gridmanager-potential/gridmanager-potential.cpp @@ -438,7 +438,7 @@ void save_text(const std::vector &cmd_units) void data_cloud(const std::vector &cmd_units) { - // data-cloud \ node|cell \ \ \ \ [\] + // data-cloud \ node|cell \ \ \ \ if (cmd_units.size() < 7) throw std::runtime_error("data-cloud: insufficient parameters."); gctl::array copy_str(7, "null"); @@ -452,33 +452,15 @@ void data_cloud(const std::vector &cmd_units) else if (copy_str[1] == "cell") d_type = ElemData; else throw std::runtime_error("open-surfer: invalid grid data type."); - text_descriptor desc; - desc.file_name_ = copy_str[2]; - desc.file_ext_ = ".txt"; - desc.col_str_ = "0,1,2"; - desc.delimiter_ = ' '; - desc.att_sym_ = '#'; - desc.head_num_ = 0; - if (copy_str[6] != "null") - { - parse_string_to_value(copy_str[6], '/', true, desc.col_str_, - desc.delimiter_, desc.att_sym_, desc.head_num_); - } + geodsv_io text_in; + text_in.load_text(copy_str[2]); - std::vector posix_vec, posiy_vec, data_vec; - read_text2vectors(desc, posix_vec, posiy_vec, data_vec); + array posi_arr; + _1d_array data_vec; + text_in.get_column_point2dc(posi_arr, 0, 1); + text_in.get_column(data_vec, 2); - array posi_arr(posix_vec.size()); - array posi_val; - - posi_val.input(data_vec); - for (size_t i = 0; i < posi_arr.size(); i++) - { - posi_arr[i].x = posix_vec[i]; - posi_arr[i].y = posiy_vec[i]; - } - - rg.load_data_cloud(posi_arr, posi_val, atof(copy_str[3].c_str()), atof(copy_str[4].c_str()), atof(copy_str[5].c_str()), copy_str[0], d_type); + rg.load_data_cloud(posi_arr, data_vec, atof(copy_str[3].c_str()), atof(copy_str[4].c_str()), atof(copy_str[5].c_str()), copy_str[0], d_type); return; } @@ -721,29 +703,9 @@ void get_profile(const std::vector &cmd_units) } else // File exist { - text_descriptor desc; - desc.file_name_ = copy_str[2]; - desc.file_ext_ = ".txt"; - desc.att_sym_ = '#'; - desc.col_str_ = "0,1"; - desc.delimiter_ = ' '; - desc.head_num_ = 0; - - if (copy_str[3] != "null") - { - gctl::parse_string_to_value(copy_str[3], '/', true, desc.col_str_, desc.delimiter_, - desc.att_sym_, desc.head_num_); - } - - std::vector xs, ys; - gctl::read_text2vectors(desc, xs, ys); - - xys.resize(xs.size()); - for (size_t i = 0; i < xys.size(); i++) - { - xys[i].x = xs[i]; - xys[i].y = ys[i]; - } + geodsv_io fio; + fio.load_text(copy_str[2]); + fio.get_column_point2dc(xys, 0, 1); rg.extract_points(copy_str[0], xys, p_data); } diff --git a/tool/gridmanager-potential/readme_gridmanager_potential.md b/tool/gridmanager-potential/readme_gridmanager_potential.md index 3256a12..4c6f266 100644 --- a/tool/gridmanager-potential/readme_gridmanager_potential.md +++ b/tool/gridmanager-potential/readme_gridmanager_potential.md @@ -22,8 +22,8 @@ Mathematic expresstions could be used to specify the grid ranges. 4. `save gmsh ` Save the grid using the Gmsh legacy format. This command will write all data that are allowed for outputing. 5. `save text ` Save the grid data to text file. -#### cloud \ node|cell \ \ \ \ [\] -Import a randomly distributed data could (namely a xyz data). Note that a grid must exist before the use this command. 'node' or 'cell' indicates whether the input could data should be node or cell registed. \ is the input data file. \,\ and \ defines an oval of which the command is used to interpolate the could data to grid points. Additional options for reading the \ could be given by the [\] argument which has a format of \/\/\/\. For example 0,1,2/,/#/0 means that the intput data could is stored int the first three columns separated by commas. The file has no head records and any line starts with '#' is an annotation. +#### cloud \ node|cell \ \ \ \ +Import a randomly distributed data could (namely a xyz data). Note that a grid must exist before the use this command. 'node' or 'cell' indicates whether the input could data should be node or cell registed. \ is the input data file. \,\ and \ defines an oval of which the command is used to interpolate the could data to grid points. #### gradient \ \ dx|dy \ Calculate directional gradient(s) of the selected data. The new gradient data will have the same