71 lines
2.1 KiB
C++
71 lines
2.1 KiB
C++
|
#include "head_func.h"
|
||
|
|
|
||
|
|
/*************************全局函数********************************/
|
||
|
|
double sign(double input){
|
||
|
|
if (input >= 0.0) return 1.0;
|
||
|
|
else return -1.0;
|
||
|
|
}
|
||
|
|
//将string转换为stringstream
|
||
|
|
stringstream str2ss(string s){
|
||
|
|
stringstream sstr;
|
||
|
|
sstr.str(""); sstr.clear(); sstr.str(s);
|
||
|
|
return sstr;
|
||
|
|
}
|
||
|
|
//测试打开输入文件 如果成功则返回0并输出信息 否则返回1
|
||
|
|
int open_infile(ifstream &infile,char* filename){
|
||
|
|
infile.open(filename);
|
||
|
|
if (!infile){
|
||
|
|
cerr << BOLDRED << "error ==> " << RESET << "file not found: " << filename << endl;
|
||
|
|
return -1;
|
||
|
|
}
|
||
|
|
return 0;
|
||
|
|
}
|
||
|
|
//测试打开输出文件 如果成功则返回0并输出信息 否则返回1
|
||
|
|
int open_outfile(ofstream &outfile,char* filename){
|
||
|
|
outfile.open(filename);
|
||
|
|
if (!outfile){
|
||
|
|
cerr << BOLDRED << "error ==> " << RESET << "fail to create the file: " << filename << endl;
|
||
|
|
return -1;
|
||
|
|
}
|
||
|
|
return 0;
|
||
|
|
}
|
||
|
|
//计算一个参考椭球或者参考球在指定纬度位置的半径
|
||
|
|
double ellipRadius(double lati,double refr,double refR){
|
||
|
|
return refr*refR/sqrt(pow(refr,2)*pow(cos(lati*Pi/180.0),2)+pow(refR,2)*pow(sin(lati*Pi/180.0),2));
|
||
|
|
}
|
||
|
|
|
||
|
|
/*************************球谐系数操作********************************/
|
||
|
|
//计算伴随勒让德函数递推系数
|
||
|
|
_2dArray get_Anm_array(int MaxOrder){
|
||
|
|
int i,j;
|
||
|
|
_2dArray cs;
|
||
|
|
cs.resize(MaxOrder);
|
||
|
|
for (i = 0; i < MaxOrder; i++)
|
||
|
|
cs[i].resize(i+1);
|
||
|
|
//向下列推计算
|
||
|
|
#pragma omp parallel for private(i,j) schedule(guided)
|
||
|
|
for (j = 0; j < MaxOrder; j++){
|
||
|
|
cs[j][j] = 0; //对角线上的值直接给0 反正用不到
|
||
|
|
for (i = j+1; i < MaxOrder; i++){
|
||
|
|
cs[i][j] = sqrt(((2.0*i-1)*(2.0*i+1))/((i-j)*(i+j)));
|
||
|
|
}
|
||
|
|
}
|
||
|
|
return cs;
|
||
|
|
}
|
||
|
|
|
||
|
|
_2dArray get_Bnm_array(int MaxOrder){
|
||
|
|
int i,j;
|
||
|
|
_2dArray cs;
|
||
|
|
cs.resize(MaxOrder);
|
||
|
|
for (i = 0; i < MaxOrder; i++)
|
||
|
|
cs[i].resize(i+1);
|
||
|
|
//向下列推计算
|
||
|
|
#pragma omp parallel for private(i,j) schedule(guided)
|
||
|
|
for (j = 0; j < MaxOrder; j++){
|
||
|
|
cs[j][j] = 0; //对角线上的值直接给0 反正用不到
|
||
|
|
for (i = j+1; i < MaxOrder; i++){
|
||
|
|
cs[i][j] = sqrt(((2.0*i+1)*(i+j-1)*(i-j-1))/((i-j)*(i+j)*(2.0*i-3)));
|
||
|
|
}
|
||
|
|
}
|
||
|
|
return cs;
|
||
|
|
}
|