diff --git a/PyGMT_out/.idea/.name b/PyGMT_out/.idea/.name new file mode 100644 index 0000000..4eb38fe --- /dev/null +++ b/PyGMT_out/.idea/.name @@ -0,0 +1 @@ +pyGMT.py \ No newline at end of file diff --git a/PyGMT_out/.idea/vcs.xml b/PyGMT_out/.idea/vcs.xml new file mode 100644 index 0000000..6c0b863 --- /dev/null +++ b/PyGMT_out/.idea/vcs.xml @@ -0,0 +1,6 @@ + + + + + + \ No newline at end of file diff --git a/PyGMT_out/pyGMT.py b/PyGMT_out/pyGMT.py index 23f8642..8525e32 100644 --- a/PyGMT_out/pyGMT.py +++ b/PyGMT_out/pyGMT.py @@ -9,7 +9,7 @@ import pygmt fig = pygmt.Figure() fig.coast(region="d",projection="Cyl_stere/12c",land="white",water="#C6E2FF",borders="1/0.2p",shorelines="1/0.2p",frame="ag") -f = open("LLA_result_for_PyGMT.txt","r") +f = open("LLA_result.txt","r") f.readline() string = f.readline() epoch = 0 diff --git a/headers/cal_mod.h b/headers/cal_mod.h index f31c251..f7e7ab9 100644 --- a/headers/cal_mod.h +++ b/headers/cal_mod.h @@ -6,17 +6,17 @@ #include"public.h" -#include"read.h" #include"2gpst.h" -#include"select_epo.h" +#include"rahcal.h" #include"xyz2blh.h" #include"degRrad.h" #include"deg2dms.h" #include"blh2enu.h" -#include"rahcal.h" +#include"type2code.h" +#include"select_epo.h" /* ------------------------------- GPS卫星位置计算函数 ------------------------------ */ -int gps_pos(int eponum, int sPRN, int best_epoch, double GPSsec, +void gps_pos(int eponum, int sPRN, int best_epoch, double GPSsec, pnav_body nav_b, pobs_head obs_h, ppos_t pos_t) { double T_S = GPSsec - pos_t[eponum].deltat[sPRN];//计算信号发射时刻 @@ -64,7 +64,7 @@ int gps_pos(int eponum, int sPRN, int best_epoch, double GPSsec, /* ------------------------------- GPS卫星位置计算&输出 ------------------------------ */ -int sat_gps_pos_clac(FILE * result_file, +void sat_gps_pos_clac(FILE * result_file, pnav_body nav_b, pobs_epoch obs_e, pobs_body obs_b, pobs_head obs_h, pstation station, ppos_t pos_t, pblh blh, penu enu, prah rah, @@ -109,7 +109,7 @@ int sat_gps_pos_clac(FILE * result_file, station[i].delta_TR == 0.0; //计算近似的信号传播时间,接收机钟差已初始化为0(伪距/光速-接收机钟差+卫星钟差) - pos_t[i].delta_t[sPRN] = (obs_b[i].obs_gps[j][C1] / C_V) - station[i].delta_TR + nav_b[best_epoch].sa0 + nav_b[best_epoch].sa1 * detat_toc + nav_b[best_epoch].sa2 * pow(detat_toc, 2); + pos_t[i].delta_t[sPRN] = (obs_b[i].obs_gps[j][Code2Type(C1C, obs_h->obstypenum_gps, obs_h->obscode_gps)] / C_V) - station[i].delta_TR + nav_b[best_epoch].sa0 + nav_b[best_epoch].sa1 * detat_toc + nav_b[best_epoch].sa2 * pow(detat_toc, 2); //判断收敛 pos_t[i].deltat[sPRN] = 0.0; @@ -171,7 +171,10 @@ int sat_gps_pos_clac(FILE * result_file, fclose(result_file); result_file = fopen(".\\Result_out\\LLA_result.txt", "a+"); fprintf(result_file, "\nG%02d[sPRN] %15.05f %15.05f %15.05f %15.05f %15.05f %15.05f %15.12f %16.13f" - ,sPRN, pos_t[i].X[sPRN], pos_t[i].Y[sPRN], pos_t[i].Z[sPRN], obs_b[i].obs_gps[j][C1], obs_b[i].obs_gps[j][10], rah->H, nav_b[best_epoch].TGD, pos_t[i].delta_clk[sPRN]); + ,sPRN, pos_t[i].X[sPRN], pos_t[i].Y[sPRN], pos_t[i].Z[sPRN] + ,obs_b[i].obs_gps[j][Code2Type(C1C, obs_h->obstypenum_gps, obs_h->obscode_gps)] + ,obs_b[i].obs_gps[j][Code2Type(C2L, obs_h->obstypenum_gps, obs_h->obscode_gps)] + ,rah->H, nav_b[best_epoch].TGD, pos_t[i].delta_clk[sPRN]); fclose(result_file); } @@ -207,7 +210,6 @@ int sat_gps_pos_clac(FILE * result_file, free(rah);rah = NULL; free(pos_t);pos_t = NULL; free(station);station = NULL; - return 1; } #endif \ No newline at end of file diff --git a/headers/public.h b/headers/public.h index 7532791..4a066d6 100644 --- a/headers/public.h +++ b/headers/public.h @@ -9,29 +9,28 @@ #define math_e 2.718281828459 #define PI 3.141592653589793 #define Earth_e 7.2921151467e-5 -#define C1 0 /* ----------------------------- GPS OBSERV TYPE ---------------------------- */ -#define C1C 101 -#define C1S 102 -#define C1L 103 -#define C1X 104 -#define C1P 105 -#define C1W 106 -#define C1Y 107 -#define C1M 108 +#define C1C 1 +#define C1S 2 +#define C1L 3 +#define C1X 4 +#define C1P 5 +#define C1W 6 +#define C1Y 7 +#define C1M 8 -#define C2C 201 -#define C2S 202 -#define C2L 203 -#define C2X 204 -#define C2P 205 -#define C2W 206 -#define C2Y 207 -#define C2M 208 +#define C2C 9 +#define C2S 10 +#define C2L 11 +#define C2X 12 +#define C2P 13 +#define C2W 14 +#define C2Y 15 +#define C2M 16 -#define C5I 501 -#define C5Q 502 -#define C5X 503 +#define C5I 17 +#define C5Q 18 +#define C5X 19 /* ------------------------------ FOR CALMOD.H ------------------------------ */ #define MAXRINEX 1000 //卫星位置结构体 @@ -70,14 +69,14 @@ double H; }blh, *pblh; //东北天结构体 - typedef struct enu + typedef struct { double E; double N; double U; }enu, *penu; //高度方位角结构体 - typedef struct rah + typedef struct { double R; double A; @@ -178,13 +177,10 @@ typedef struct obs_head double ANTH;//天线高 double ANTdeltaE;//天线中心对于测站标志在东方向上的偏移量 double ANTdeltaN;//天线中心对于测站标志在北方向上的偏移量 - int WAFL1;//缺省的L1载波的波长因子 - int WAFL2;//缺省的L2载波的波长因子 - int WAFflag; int obstypenum_gps;//GPS观测值类型数量 int obstypenum_bds; - char obstype_gps[36];//rinex3.04GPS观测值类型列表 - char obstype_bds[36]; + int obscode_gps[36];//rinex3.04GPS观测值类型列表 + int obscode_bds[36]; double interval;//观测值的历元间隔 int f_y;//第一个观测记录的时刻 int f_m; diff --git a/headers/read.h b/headers/read.h index 65cdada..2728264 100644 --- a/headers/read.h +++ b/headers/read.h @@ -1,10 +1,10 @@ #ifndef _READ_H_ #define _READ_H_ - #include #include #include"public.h" +#include"type2code.h" /* ------------------------ 将字符串转换为浮点数,i起始位置,n输入多少个字符 ----------------------- */ static double strtonum(const char* buff, int i, int n) { @@ -84,13 +84,6 @@ extern void read_o_h(FILE* fp_obs, pobs_head obs_h) obs_h->ANTdeltaN = strtonum(buff, 14 + 14, 14); continue; } - else if (strstr(lable, "WAVELENGTH FACT L1/2")) - { - obs_h->WAFL1 = (int)strtonum(buff, 0, 6); - obs_h->WAFL2 = (int)strtonum(buff, 6, 6); - obs_h->WAFflag = (int)strtonum(buff, 6 + 6, 6); - continue; - } else if (strstr(lable, "SYS / # / OBS TYPES")) { char flag = {0}; @@ -102,19 +95,19 @@ extern void read_o_h(FILE* fp_obs, pobs_head obs_h) { for (i = 0; i < obs_h->obstypenum_gps; i++) { - strncpy(&(obs_h->obstype_gps[i]), buff + 7 + 4 * i, 3); + obs_h->obscode_gps[i] = Type2Code(i, buff); } } else if (obs_h->obstypenum_gps > 13) { for (i = 0; i < 13; i++) { - strncpy(&(obs_h->obstype_gps[i]), buff + 7 + 4 * i, 3); + obs_h->obscode_gps[i] = Type2Code(i, buff); } fgets(buff, MAXRINEX, fp_obs); for (i = 0; i < obs_h->obstypenum_gps - 13; i++) { - strncpy(&(obs_h->obstype_gps[i + 13]), buff + 7 + 4 * i, 3); + obs_h->obscode_gps[i + 13] = Type2Code(i, buff); } } continue; @@ -126,19 +119,19 @@ extern void read_o_h(FILE* fp_obs, pobs_head obs_h) { for (i = 0; i < obs_h->obstypenum_bds; i++) { - strncpy(&(obs_h->obstype_bds[i]), buff + 7 + 4 * i, 3); + obs_h->obscode_bds[i] = Type2Code(i, buff); } } else if (obs_h->obstypenum_bds > 13) { for (i = 0; i < 13; i++) { - strncpy(&(obs_h->obstype_bds[i]), buff + 7 + 4 * i, 3); + obs_h->obscode_bds[i] = Type2Code(i, buff); } fgets(buff, MAXRINEX, fp_obs); for (i = 0; i < obs_h->obstypenum_bds - 13; i++) { - strncpy(&(obs_h->obstype_bds[i + 13]), buff + 7 + 4 * i, 3); + obs_h->obscode_bds[i + 13] = Type2Code(i, buff); } } continue; @@ -213,7 +206,7 @@ extern void read_o_b(FILE* fp_obs, pobs_epoch obs_e, pobs_body obs_b, int type_g else if (flag != 'G')//获取flag标识,判断GPS卫星为否 { //判断是否为BDS卫星 - if (flag == 'C')//获取flag标识,判断BDS卫星为真,读取该卫星PRN号 + if (flag == 'C')//获取flag标识,判断BDS卫星为真,读取该卫星PRN号 { obs_e[n].sPRN[i] = strtonum(buff, 1, 2);//将卫星编号存储到第n个历元中,第i个PRN位置 obs_e[n].sPRN_BDS[j_c] = obs_e[n].sPRN[i];//将BDS卫星编号与总观测卫星号关联 @@ -232,11 +225,9 @@ extern void read_o_b(FILE* fp_obs, pobs_epoch obs_e, pobs_body obs_b, int type_g strncpy(&flag, buff + 0, 1); } } - } n++; - } - + } } /* -------------------------------------------------------------------------- */ @@ -244,7 +235,7 @@ extern void read_o_b(FILE* fp_obs, pobs_epoch obs_e, pobs_body obs_b, int type_g /* -------------------------------------------------------------------------- */ /* ---------------------- 获取文件数据块行数,从END OF HEADER后开始起算 --------------------- */ -extern int getgpssatnum(FILE* fp_nav) +int getgpssatnum(FILE* fp_nav) { int gps_satnum = 0; int flag = 0; @@ -269,8 +260,7 @@ extern int getgpssatnum(FILE* fp_nav) gps_satnum++; break; } - } - + } } if (strstr(lable, "END OF HEADER")) diff --git a/headers/type2code.h b/headers/type2code.h new file mode 100644 index 0000000..5864813 --- /dev/null +++ b/headers/type2code.h @@ -0,0 +1,59 @@ +#ifndef _TYPE2CODE_H_ +#define _TYPE2CODE_H_ + +#include +#include + +#include"public.h" + +/* -------------------------------------------------------------------------- */ +/// @brief 将观测值类型转换为PUBLIC.H中定义的数字代码 +/// @param i 读取循环中传入的循环次数 +/// @param buff 传入的观测值类型 +/// @return 返回观测值代码 +/* -------------------------------------------------------------------------- */ +int Type2Code(int i, char buff[1000]) +{ + int code_return = 0; + char code_contrast[4]; + char code_original[4]; + + for(int j = 0; j < 19; j++) + { + memset(code_contrast,0,sizeof(code_contrast)); + memset(code_original,0,sizeof(code_original)); + const char codeindex[]="C1CC1SC1LC1XC1PC1WC1YC1M" + "C2CC2SC2LC2XC2PC2WC2YC2M" + "C5IC5QC5X"; + strncpy(code_contrast, codeindex + 0 + j * 3 , 3); + strncpy(code_original, buff + 7 + 4 * i , 3); + if(strcmp(code_contrast,code_original) == 0) + { + code_return = j + 1; + break; + } + } + return code_return; +} + +/* -------------------------------------------------------------------------- */ +/// @brief 将读取的观测值类型匹配类型代码存储数组 +/// @param typecode 希望匹配的观测值 +/// @param typenum 观测值类型总数量 +/// @param typearr 类型代码存储数组 +/// @return 返回对应的观测值类型代码的数组位置 +/* -------------------------------------------------------------------------- */ +int Code2Type(int typecode, int typenum, int typearr[36]) +{ + int pos; + for(pos = 0; pos < typenum; pos++) + { + if(typearr[pos] == typecode) + { + break; + } + } + return pos; +} + +#endif \ No newline at end of file