logplus/BaseFun/src/BaseFun.cpp
2025-10-29 17:23:30 +08:00

3157 lines
71 KiB
C++
Raw Blame History

This file contains ambiguous Unicode characters

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

// BaseFun.cpp : 定义 DLL 的初始化例程。
//
#include "BaseFun.h"
#include "math.h"
#include <QDir>
//#include "ObjectArchive.h"
#include <QCoreApplication>
#include <QTextStream>
#include <QTextCodec>
#include <QLibrary>
#define REPR_INT 1
#define REPR_SHORT 2
#define REPR_LONG 3
#define REPR_FLOAT 4
#define REPR_DOUBLE 5
#define REPR_STRING 6
#define REPR_CHAR 7
#define REPR_UCHAR 8
#define REPR_USHORT 9
#define REPR_UINT 10
#define REPR_ULONG 11
void removeStr(char* sdes,QString b,QString e)
{
QString des=sdes;
removeStr(des,b,e);
strcpy(sdes,des.toStdString().c_str());
}
void removeStr(QString &des,QString b,QString e)
{
while(1){
int ind1=des.indexOf(b);
if(ind1>-1) {
int indx=des.indexOf(e);
if(indx<0) des=des.left(ind1);
else {
des.remove(ind1,indx-ind1+1);
}
}
else break;
}
}
void CalTextAlign(int & hAlign, int &vAlign,bool isRot)
{
switch(hAlign)
{
case 1:
if(!isRot)
hAlign = Qt::AlignLeft;
else hAlign = Qt::AlignBottom;
break;
case 2:
if(!isRot)
hAlign = Qt::AlignRight;
else hAlign = Qt::AlignTop;
break;
case 0:
default:
if(!isRot)
hAlign = Qt::AlignHCenter;
else hAlign = Qt::AlignVCenter;
break;
}
switch(vAlign)
{
case 0:
if(!isRot)
vAlign = Qt::AlignTop;
else vAlign = Qt::AlignLeft;
break;
case 2:
if(!isRot)
vAlign = Qt::AlignBottom;
else vAlign = Qt::AlignRight;
break;
case 1:
default:
if(!isRot)
vAlign = Qt::AlignVCenter;
else vAlign = Qt::AlignHCenter;
break;
}
}
void CalTextWrap( QString &SrcText,double width,QFont font,QString &OutText,double &height,bool IsWellSectonHorizonLayout )
{
int len=SrcText.length();
if (len==0)
return;
QFontMetrics mt(font);
QString temp,temps;
QStringList tem;
for(int i=0;i<len;i++)
{
temps=temp;
temp.append(SrcText.at(i));
float w=mt.width(temp);
float h=mt.height();
if(w>width) {
temp=SrcText.at(i);
tem.append(temps);
}
if(i==len-1) {
tem.append(temp);
}
}
OutText=tem.join(" \n");
}
#include <iostream>
#include <string>
#include <stack>
#include <cmath>
#include <stdexcept>
#include <cctype>
using namespace std;
/*
class FormulaParser {
private:
string expression;
size_t index;
// 跳过空白字符
void skipWhitespace() {
while (index < expression.size() && isspace(expression[index])) {
index++;
}
}
// 获取当前字符
char peek() {
skipWhitespace();
return index < expression.size() ? expression[index] : '\0';
}
// 消费当前字符并前进
char consume() {
skipWhitespace();
return index < expression.size() ? expression[index++] : '\0';
}
// 解析数字
double parseNumber() {
string numStr;
while (index < expression.size() &&
(isdigit(expression[index]) || expression[index] == '.')) {
numStr += expression[index++];
}
return stod(numStr);
}
// 解析因子(数字或括号表达式)
double parseFactor()
{
char next = peek();
if (next == '(') {
consume(); // 消费 '('
double result = parseExpression();
if (peek() != ')') {
throw runtime_error("Missing closing parenthesis");
}
consume(); // 消费 ')'
return result;
} else if (isdigit(next)) {
return parseNumber();
} else if (next == '+' || next == '-') {
consume();
double factor = parseFactor();
return next == '+' ? factor : -factor;
} else {
throw runtime_error("Unexpected character in expression");
}
}
// 解析指数运算
double parseExponent() {
double left = parseFactor();
while (true) {
char op = peek();
if (op == '^') {
consume();
double right = parseFactor();
left = pow(left, right);
} else {
break;
}
}
return left;
}
// 解析乘除运算
double parseTerm() {
double left = parseExponent();
while (true) {
char op = peek();
if (op == '*' || op == '/') {
consume();
double right = parseExponent();
if (op == '*') {
left *= right;
} else {
if (right == 0) {
left=0;
// throw runtime_error("Division by zero");
}
else left /= right;
}
} else {
break;
}
}
return left;
}
// 解析加减运算
double parseExpression() {
double left = parseTerm();
while (true) {
char op = peek();
if (op == '+' || op == '-') {
consume();
double right = parseTerm();
if (op == '+') {
left += right;
} else {
left -= right;
}
} else {
break;
}
}
return left;
}
public:
FormulaParser(const string& expr) : expression(expr), index(0) {}
double evaluate() {
double result = parseExpression();
if (peek() != '\0') {
throw runtime_error("Unexpected characters at end of expression");
}
return result;
}
};
int main1() {
while (true) {
cout << "Enter an expression (or 'quit' to exit): ";
string input;
getline(cin, input);
if (input == "quit") {
break;
}
try {
FormulaParser parser(input);
double result = parser.evaluate();
cout << "Result: " << result << endl;
} catch (const exception& e) {
cerr << "Error: " << e.what() << endl;
}
}
return 0;
}
*/
#include <iostream>
#include <string>
#include <vector>
#include <map>
#include <cmath>
#include <stdexcept>
#include <cctype>
#include <functional> // 添加function头文件
using namespace std;
class FormulaParser {
public:
FormulaParser(const string& expression, const map<string, double>& variables)
: expr(expression), pos(0), variables(variables) {
// 预定义函数
functions["sin"] = [](double x) { return sin(x); };
functions["cos"] = [](double x) { return cos(x); };
functions["tan"] = [](double x) { return tan(x); };
functions["log"] = [](double x) { return log(x); };
functions["exp"] = [](double x) { return exp(x); };
functions["sqrt"] = [](double x) { return sqrt(x); };
// 可以添加更多函数...
functions["abs"] = [](double x) { return abs(x); };
}
double parse() {
pos = 0;
double result = parseExpression();
if (pos < expr.size()) {
throw runtime_error("Unexpected character at position " + to_string(pos));
}
return result;
}
private:
string expr;
size_t pos;
map<string, double> variables;
map<string, function<double(double)>> functions;
void skipWhitespace() {
while (pos < expr.size() && isspace(expr[pos])) {
pos++;
}
}
char peek() {
skipWhitespace();
return pos < expr.size() ? expr[pos] : '\0';
}
char consume() {
skipWhitespace();
return pos < expr.size() ? expr[pos++] : '\0';
}
double parseExpression() {
double left = parseTerm();
while (true) {
char op = peek();
if (op == '+' || op == '-') {
consume();
double right = parseTerm();
if (op == '+') {
left += right;
} else {
left -= right;
}
} else {
break;
}
}
return left;
}
double parseTerm() {
double left = parseFactor();
while (true) {
char op = peek();
if (op == '*' || op == '/') {
consume();
double right = parseFactor();
if (op == '*') {
left *= right;
} else {
if (right == 0) {
throw runtime_error("Division by zero");
}
left /= right;
}
} else {
break;
}
}
return left;
}
double parseFactor() {
double left = parsePrimary();
while (peek() == '^') {
consume();
double right = parsePrimary();
left = pow(left, right);
}
return left;
}
double parsePrimary() {
char c = peek();
if (c == '(') {
consume();
double result = parseExpression();
if (peek() != ')') {
throw runtime_error("Expected ')'");
}
consume();
return result;
} else if (isalpha(c)) {
return parseIdentifier();
} else if (isdigit(c) || c == '.') {
return parseNumber();
} else if (c == '-') {
consume();
return -parsePrimary();
} else {
throw runtime_error("Unexpected character: " + string(1, c));
}
}
double parseNumber() {
string numStr;
while (isdigit(peek())) {
numStr += consume();
}
if (peek() == '.') {
numStr += consume();
while (isdigit(peek())) {
numStr += consume();
}
}
// 科学计数法支持
if (peek() == 'e' || peek() == 'E') {
numStr += consume();
if (peek() == '+' || peek() == '-') {
numStr += consume();
}
while (isdigit(peek())) {
numStr += consume();
}
}
return stod(numStr);
}
double parseIdentifier() {
string id;
while (isalnum(peek())) {
id += consume();
}
// 检查是否是函数调用
if (peek() == '(') {
consume();
double arg = parseExpression();
if (peek() != ')') {
throw runtime_error("Expected ')' after function argument");
}
consume();
auto it = functions.find(id);
if (it != functions.end()) {
return it->second(arg);
} else {
throw runtime_error("Unknown function: " + id);
}
} else {
// 变量处理
auto it = variables.find(id);
if (it != variables.end()) {
return it->second;
} else {
throw runtime_error("Unknown variable: " + id);
}
}
}
};
float lst(int code,int point,float *data,char *fieldname,char *formu,float *value,float rlev)
{
int i,k=0,in1,in2,in3;
float br1,br2,br3;
if(point<=0) point=1;
switch(code) {
case 1:
br1=9999.;
for(i=0;i<point;i++)
{
if(data[i+k]<br1) br1=data[i+k];
}
return br1;
break;
case 2:
for(br1=0.0,i=0;i<point;i++) br1+=data[k+i];
br1=br1/point;
for(in2=0,br2=0.,i=0;i<point;i++) if(data[i+k]<=br1) {
in2++;
br2+=data[k+i];
}
if (in2 <= 0)return br1;
br2=br2/in2;
for(in3=0,br3=0,i=0;i<point;i++) if(data[i+k]<=br2) {
in3++;
br3+=data[i+k];
}
if (in3 <= 0)return br2;
br3=br3/in3;
if (in3>=(float)point/3. )return br3;
else return br2;
break;
case 3:
for(br1=0,i=0;i<point;i++) br1+=data[i+k];
br1/=point;
return br1;
case 4:
for(br1=0.,i=0;i<point;i++) br1+=data[i+k];
br1/=point;
for(in2=0,br2=0.,i=0;i<point;i++) if(data[k+i]>=br1) {
in2++; br2+=data[k+i];
}
if(in2 <= 0.0) return br1;
br2/=in2;
for(in3=0,br3=0.,i=0;i<point;i++) if(data[i+k]>=br2) {
in3++;
br3+=data[k+i];
}
if (in3 <=0.0) return br2;
br3/=in3;
if ( in3>= (float)point/3. ) return br3;
else return br2;
break;
case 5:
br1=-9999.;
for(i=0;i<point;i++) if(data[i+k]>br1) br1=data[i];
return br1;
break;
case 6:
br1=-9999.;
if(point>0) br1=data[0];
return br1;
break;
case 7:
br1=-9999.;
if(point>0) br1=data[point-1];
return br1;
break;
case 8:
{
br1=-9999.;
QString newFormu=formu;
newFormu=newFormu.toUpper();
QString Temp=fieldname;
QStringList Field,sss;
sss =Temp.split(",");
for(int i=0;i<sss.size();i++){
Field.append(sss[i].split(" "));
}
Field.removeAll("");
/*
qsort(Field,data,0,Field.size()-1,1);
char buf[200];
for(int i=0;i<Field.size();i++)
{
int pos=newFormu.indexOf(Field[i],0,Qt::CaseInsensitive);
if(pos>-1) {
// if(pos>0&&newFormu[pos-1].isLetter()) continue;
// if(pos+Field.length()<newFormu.length()&&newFormu[pos+Field.length()].isLetter()) continue;
}
else continue;
sprintf(buf,"%f",data[i]);
newFormu.replace(Field[i].toUpper(),buf);
}
*/
map<string, double> vars;
for(int i=0;i<Field.size();i++)
{
vars[Field[i].toUpper().toStdString()]=data[i];
}
try
{
FormulaParser parser(newFormu.toStdString(),vars);
br1 = parser.parse();
}
catch (const exception& e) {
cerr << "Error: " << e.what() << endl;
}
return br1;
}
break;
case 9:
br1=-9999.;
if(point>0) br1=data[1];
return br1;
break;
case 10:
br1=-9999.;
{
int first=0;
int end=point;
if(point*rlev>2) {
first=2;
end-=2;
}
float val=0;
for(int i=first;i<end;i++) {
val+=data[i];
}
int count=end-first;
if(count>0) br1=val/count;
}
return br1;
break;
}
return -9999;
}
char m_strValue[200];
QStringList GetStringList(QString line,bool IsSpa,bool IsTab,bool IsCom,bool IsSem,bool DelDubSpa)
{
if(!IsCom)line.replace(",","");//如果逗号不是分隔符,清除逗号
if(!IsSpa)line.replace(" ","");
else
{
//连续空格视为单个处理
if(DelDubSpa)while(line.indexOf(" ")>=0)line.replace(" "," ");
line=line.simplified();
//逗号空格连续存在的只留一个??
while(line.indexOf(" ,")>=0)line.replace(" ,",",");//line.replace(" ,"," ");2020.5.9
while(line.indexOf(", ")>=0)line.replace(", ",",");//line.replace(", "," ");
line.replace(" ",",");
}
if(!IsTab)line.replace(" ","");//tab
else line.replace(" ",",");
if(!IsSem)line.replace(";","");
else line.replace(";",",");
QStringList List=line.split(",");
return List;
}
double GetData(int repCode,char *buffer,int repLen)
{
double yy=-99999;
if(!buffer)
{
return 0;
}
switch(repCode)
{
case REPR_INT: //0
yy=(double)(*((int*)buffer));
break;
case REPR_SHORT: //1
yy=(double)(*((short *)buffer));
break;
case REPR_LONG://2
yy=(double)(*((long *)buffer));
break;
case REPR_FLOAT://3
yy=(double)(*((float *)buffer));
break;
case REPR_DOUBLE://4
yy=(double)(*((double *)buffer));
break;
case REPR_CHAR://5
yy=(double)(*((char *)buffer));
break;
case REPR_UCHAR://6
yy=(double)(*((unsigned char *)buffer));
break;
case REPR_USHORT://7
yy=(double)(*((unsigned short *)buffer));
break;
case REPR_UINT://8
yy=(double)(*((unsigned int *)buffer));
break;
case REPR_ULONG://9
yy=(double)(*((unsigned long *)buffer));
break;
case REPR_STRING://10
if(repLen>=200) repLen=199;
if(repLen>0) {
memmove(m_strValue,buffer,repLen);
m_strValue[repLen]=0;
}
else m_strValue[0]=0;
break;
}
return yy;
}
void SetData(int repCode,char *buffer,double yy)
{
if(!buffer)
{
return;
}
switch(repCode)
{
case REPR_INT: //0
(*((int*)buffer))=(int)yy;
break;
case REPR_SHORT: //1
(*((short *)buffer))=(short)yy;
break;
case REPR_LONG://2
(*((long *)buffer))=(long)yy;
break;
case REPR_FLOAT://3
(*((float *)buffer))=yy;
break;
case REPR_DOUBLE://4
(*((double *)buffer))=(double)yy;
break;
case REPR_CHAR://5
(*((char *)buffer))=(char)yy;
break;
case REPR_UCHAR://6
(*((unsigned char *)buffer))=(unsigned char)yy;
break;
case REPR_USHORT://7
(*((unsigned short *)buffer))=(unsigned short)yy;
break;
case REPR_UINT://8
(*((unsigned int *)buffer))=(unsigned int)yy;
break;
case REPR_ULONG://9
(*((unsigned long *)buffer))=(unsigned long )yy;
break;
case REPR_STRING://10
// *yy=-99999;
break;
}
return;
}
int akima (float *x,float *y,int nx,float *c)
{
int i, no;
float t1, t2, b, rm1, rm2, rm3, rm4;
if (nx < 4) return (-1);
for (i = 1; i < nx; i++)
{
if (x[i] <= x[i-1])
{
return (-2);
}
}
rm3 = (y[1] - y[0])/(x[1] - x[0]);
t1 = rm3 - (y[1] - y[2])/(x[1] - x[2]);
rm2 = rm3 + t1;
rm1 = rm2 + t1;
/* get slopes */
no = nx - 2;
for (i = 0; i < nx; i++) {
if (i >= no)
rm4 = rm3 - rm2 + rm3;
else
rm4 = (y[i+2] - y[i+1])/(x[i+2] - x[i+1]);
t1 = (float)fabs(rm4 - rm3);
t2 = (float)fabs(rm2 - rm1);
b = t1 + t2;
c[3*i] =(float) ((b != 0.0) ? (t1*rm2 + t2*rm3) / b : 0.5*(rm2 + rm3));
rm1 = rm2;
rm2 = rm3;
rm3 = rm4;
}
no = nx - 1;
/* compute the coefficients for the nx-1 intervals */
for (i = 0; i < no; i++) {
t1 =(float) (1.0 / (x[i+1] - x[i]));
t2 = (y[i+1] - y[i])*t1;
b = (c[3*i] + c[3*i+3] - t2 - t2)*t1;
c[3*i+2] = b*t1;
c[3*i+1] = -b + (t2 - c[3*i])*t1;
}
return (0);
}
/*---------------------------------------------------------------------------------------------------------------
*/
/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
* CUBSPL computes the coefficients for a natural cubic spline. The
* algorithm is the same as the one in the IMSL library
* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
*/
int cubspl (float x[],float y[],int nx,float c[])
{
int im1, i, jj, j, mm1, mp1, m, nm1, nm2;
float divdf1, divdf3, dtau, g, cnx[3];
if (nx < 4) return (-1);
nm1 = nx - 1;
for (m = 1; m < nm1; m++) {
mm1 = m - 1;
c[3*m+1] = x[m] - x[mm1];
if (c[3*m+1] <= 0.0) return (-2);
c[3*m+2] = (y[m] - y[mm1])/c[3*m+1];
}
cnx[1] = x[nx-1] - x[nm1-1];
if (cnx[1] <= 0.0) return (-3);
cnx[2] = (y[nx-1] - y[nm1-1])/cnx[1];
nm2 = nx - 2;
c[2] = c[7];
c[1] = c[4] + c[7];
c[0] =(float) (((c[4] + 2.*c[1])*c[5]*c[7] + c[4]*c[4]*c[8])/c[1]);
for (m = 1; m < nm2; m++) {
mp1 = m + 1;
mm1 = m - 1;
g = -c[3*mp1+1]/c[3*mm1+2];
c[3*m] =(float) (g*c[3*mm1] + 3.*c[3*m+1]*c[3*mp1+2] + 3.*c[3*mp1+1]*c[3*m+2]);
c[3*m+2] =(float) (g*c[3*mm1+1] + 2.*c[3*m+1] + 2.*c[3*mp1+1]);
}
g = -cnx[1]/c[3*nm2-1];
c[3*nm1-3] =(float) (g*c[3*nm2-3] + 3.*c[3*nm1-2]*cnx[2] + 3.*cnx[1]*c[3*nm1-1]);
c[3*nm1-1] =(float) (g*c[3*nm2-2] + 2.*c[3*nm1-2] + 2.*cnx[1]);
g = c[3*nm1-2] + cnx[1];
cnx[0] =(float) (((cnx[1] + 2.*g)*cnx[2]*c[3*nm1-2] + cnx[1]*cnx[1]*(y[nm1-1]-y[nx-3])/c[3*nm1-2])/g);
g = -g/c[3*nm1-1];
cnx[2] = c[3*nm1-2] + g*c[3*nm1-2];
cnx[0] = (g*c[3*nm1-3] + cnx[0])/cnx[2];
c[3*nm1-3] = (c[3*nm1-3]-c[3*nm1-2]*cnx[0])/c[3*nm1-1];
for (jj = 0; jj < nm2; jj++) {
j = nm1 - jj - 1;
c[3*j] = (c[3*j]-c[3*j+1]*c[3*j+3])/c[3*j+2];
}
for (i = 1; i < nm1; i++) {
im1 = i-1;
dtau = c[3*i+1];
divdf1 = (y[i]-y[im1])/dtau;
divdf3 =(float) (c[3*im1] + c[3*i] - 2.*divdf1);
c[3*im1+1] = (divdf1 - c[3*im1] - divdf3)/dtau;
c[3*im1+2] = divdf3/(dtau*dtau);
}
dtau = cnx[1];
divdf1 = (y[nx-1]-y[nm1-1])/dtau;
divdf3 =(float) (c[3*nm1-3] + cnx[0] - 2.*divdf1);
c[3*nm1-2] = (divdf1-c[3*nm1-3]-divdf3)/dtau;
c[3*nm1-1] = divdf3/(dtau*dtau);
return (0);
}
/****************************************************************/
/** this program was refered to "Numerical Recipes in C
The Art of Scientific Computing Second Edition "
Wroted by W. H. Press, S. A. Teukolsky, W. T. Vetterling,
B. P. Flannery
**/
/**************************************************************/
int spline(float x[],float y[], int n,float yp1, float ypn,float y2[])
/* x[],y[]: input point data;
n: data bumber of input data;
yp1: the first rank differential of point 1;
ypn: the first rank differential of point 2;
*/
{
int i,k;
float p,qn,sig,un,*u;
if(n<4) return(-1);
u=(float *)malloc(sizeof(float)*n);
if(yp1>0.99e30) y2[0]=u[0]=0.0; /** the low boundery condition set to zero */
else {
y2[0]=-0.5;
u[0]=(float)((3.0/(x[1]-x[0]))*((y[1]-y[0])/(x[1]-x[0])-yp1));
}
for(i=1;i<n-1;i++){
sig=(x[i]-x[i-1])/(x[i+1]-x[i-1]);
p=sig*y2[i-1]+(float)2.0;
y2[i]=(sig-(float)1.0)/p;
u[i]=(y[i+1]-y[i])/(x[i+1]-x[i])-(y[i]-y[i-1])/(x[i]-x[i-1]);
u[i]=((float)6.0*u[i]/(x[i+1]-x[i-1])-sig*u[i-1])/p;
}
if(ypn>0.99e30) qn=un=0.0;
else {
qn=0.5;
un=((float)3.0/(x[n-1]-x[n-2]))*(ypn-(y[n-1]-y[n-2])/(x[n-1]-x[n-2]));
}
y2[n-1]=(un-qn*u[n-2])/(qn*y2[n-2]+(float)1.0);
for(k=n-2;k>=0;k--)
y2[k]=y2[k]*y2[k+1]+u[k];
free((char *)u);
return(0);
}
char *memory (char *prev_addr, int n, int size, char *progname)
{
char *tmp;
if (prev_addr) {
if ((tmp = (char *)realloc ((char *) prev_addr, (unsigned) (n * size))) == NULL) {
fprintf (stderr, "%s: Could not reallocate more memory, n = %d\n", progname, n);
return 0;
}
}
else {
if ((tmp = (char *)(char *)calloc ((unsigned) n, (unsigned) size)) == NULL) {
fprintf (stderr, "%s: Could not allocate memory, n = %d\n", progname, n);
return 0;
}
}
return (tmp);
}
int intpol (float *x,float *y,int n,int m,float *u,float *v,int mode)
{
int i, j, err_flag = 0;
float dx, *c, nanq, zip = 0.0;
float h,a,b;
char *CNULL = (char *)NULL;
if (n < 4) mode = 0;
char strTmp[7] = "intpol";
if (mode > 0) c = (float *) memory (CNULL, 3*n, sizeof(float), strTmp);
if (mode == 0) { /* Linear interpolation */
for (i = 1; i < n; i++)
if (x[i] <= x[i-1])
return (-1);
}
else if (mode == 1) /* Akimas spline */
err_flag = akima(x,y,n,c);
else if (mode == 2) /* Natural cubic spline */
/* err_flag = cubspl(xx,yy,n,c);
*/
err_flag=spline(x,y, n,(float)0.99e34,(float)0.99e34,c);
else if(mode==3) err_flag = cubspl(x,y,n,c);
else
err_flag = -4;
if (err_flag != 0) return (err_flag);
/* Compute the interpolated values from the coefficients */
j = 0;
nanq =(float) 1.0e-30;
for (i = 0; i < m; i++) {
if (u[i] < x[0] || u[i] > x[n-1]) { /* Desired point outside data range */
v[i] = nanq;
continue;
}
while (x[j] <= u[i] && j < n) j++;
if (j > 0) j--;
dx = u[i] - x[j];
if (mode == 0)
v[i] = (y[j+1]-y[j])*dx/(x[j+1]-x[j]) + y[j];
else if(mode ==2){
h=x[j+1]-x[j];
if(h==0.0) {fprintf(stderr,"Bad xa input ro routing splint\n");
v[i]=y[j];
}
else {a=(x[j+1]-u[i])/h;
b=(u[i]-x[j])/h;
v[i]=a*y[j]+b*y[j+1]+((a*a*a-a)*c[j]+(b*b*b-b)*c[j+1])*(h*h)/(float)6.0;
/*printf("a=%f b=%f v[%d]=%f\n",a,b,i,v[i]); */
}
}
else
v[i] = ((c[3*j+2]*dx + c[3*j+1])*dx + c[3*j])*dx + y[j];
}
if (mode > 0) free ((char *)c);
return (0);
}
///360度插值
int Resamples(char *Data0,char *Data1,int N0,int N1,int RepCode,int CodeLen)
{
float *InData=new float[N0];
float *OutData=new float[N1];
for(int i=0;i<N0;i++) {
InData[i]=GetData(RepCode,&Data0[i*CodeLen],CodeLen);
}
float ino=N0;
float ono=N1;
//Resamples(InData,OutData,ino,ono);
int pos=0;
float val1=0,val2=0,val=0;
for(int i=0;i<N1;i++) {
pos=i*N0/N1;
// val1=InData[pos];
// val2=InData[ps+1];
// Data1[i]=(i-i*N0/N1)/(1)*(val2-val1)+val1;
OutData[i]=InData[pos];
}
for(int i=0;i<N1;i++) {
SetData(RepCode,&Data1[i*CodeLen],OutData[i]);
}
delete InData;
delete OutData;
return 1;
}
int Resamples(float *Data0,float *Data1,float N0,float N1)
{
//Data0是原来的N0道数据
//Data1为插好的N1道数据
float *data0= (float*)malloc((N0+1)*sizeof(float));
float *x0=(float*)malloc((N0+1)*sizeof(float));
float *x1=(float*)malloc(N1*sizeof(float));
for(int i=0;i<N0+1;i++)
{
x0[i]=(i*1.0)*360.0/(1.0*N0);
if(i<N0)
data0[i]=Data0[i];
else
data0[i]=Data0[0];
}
for(int i=0;i<N1;i++)
{
x1[i]=(i*1.0)*360.0/(1.0*N1);
}
intpol(x0,data0,N0+1,N1,x1,Data1,1);
free(data0);
free(x0);
free(x1);
return 1;
}
void LeastSquare(int size,float *x, float* y,float &a,float &b)
{
double t1=0, t2=0, t3=0, t4=0;
for(int i=0; i<size; ++i)
{
t1 += x[i]*x[i];
t2 += x[i];
t3 += x[i]*y[i];
t4 += y[i];
}
a = (t3*size - t2*t4) / (t1*size - t2*t2);
//b = (t4 - a*t2) / x.size();
b = (t1*t4 - t2*t3) / (t1*size - t2*t2);
}
float GetMaxFreq(float m_min,float m_max,float *inFreq,int count)
{
float x=0;
int j=-1;
for(int i=0;i<count;i++) {
if(inFreq[i]>x) {
x=inFreq[i];
j=i;
}
}
float val=0;
if(j!=-1) val=j*(m_max-m_min)/count;
return val;
}
int GetFreq(float *pVal,int dataSize,float m_min,float m_max,int count,float *outFreq)
{
for (int j = 0; j < count; ++j)
{
outFreq[j] = 0;
}
//进行统计
int effectiveSize = 0;
//等于最大值得放入最后一个区间
for (int i = 0; i < (int)dataSize; ++i)
{
if (qAbs(pVal[i] - m_max) < 1e-10)
{
outFreq[count - 1]+= 1;
effectiveSize += 1;
continue;
}
int ss=(pVal[i]-m_min)/((m_max - m_min)/(float)count);
if(ss<0) continue;
if(ss>=count) continue;
outFreq[ss] += 1;
effectiveSize += 1;
}
return effectiveSize;
}
float mean(vector<float>resultSet)
{
std::vector<float>::iterator d=resultSet.begin();
double sum=0;
for (;d!=resultSet.end();d++) {
sum+=*d;
}
float mean1 = sum / resultSet.size(); //均值
return mean1;
}
float stdev(vector<float>resultSet)
{
double mean1=mean(resultSet);
double accum = 0.0;
std::vector<float>::iterator d=resultSet.begin();
for (;d!=resultSet.end();d++)
{
accum += (*d-mean1)*(*d-mean1);
};
double stdev1 = sqrt(accum/(resultSet.size()-1)); //方差
return stdev1;
}
int IsNumberic(QString src)
{
QByteArray ba = src.toLatin1();
const char *s = ba.data();
while(*s)
{
if(*s>='0'&& *s<='9')
s++;
else if(*s == '.')
s++;
else if(*s == '-')
s++;
else
break;
}
if(*s)
return -1;//不是纯数字
else
return 0;//纯数字
}
int binarySearch(float*darray, int n, float data,float err)//
{
int low, mid, high;
if(darray == NULL) return -1;//
low = 0;
high = n - 1;
while(low <= high) {
mid = (low + high) / 2;
if(fabs(data-darray[mid])<err) return mid; //
else if(data < darray[mid]) high = mid - 1;
else low = mid + 1;
}
return -1;//
}
int bSearch(double*darray, int n, double data,double err)//
{
int low, mid, high;
if(darray == NULL) return -1;//
low = 0;
high = n - 1;
while(low <= high) {
mid = (low + high) / 2;
if(fabs(data-darray[mid])<err) return mid;
if(mid+1<n&&data>darray[mid]&&data<darray[mid+1]) return mid; //
else if(data < darray[mid]) high = mid - 1;
else low = mid + 1;
}
return -1;//
}
//bool FLOATPROPERTY::Serialize( CObjectArchive &ar )
// {
// if(ar.IsStoring())
// {
// BEGIN_WRITE_OBJECT( ar , 1 );
// BEGIN_WRITE_BLOCK( ar , 1 );
// ar << m_size;
// ar.writeRawData((char *)m_vProperty,m_size*sizeof(float));
// END_WRITE_BLOCK( ar , 1);
// END_WRITE_OBJECT( ar );
// }
// else
// {
// BEGIN_READ_OBJECT( ar,1 );
// BEGIN_READ_BLOCK( 1 );
// ar >> m_size;
// if(m_vProperty) delete m_vProperty;
// m_vProperty =new float[m_size];//(float *)malloc(sizeof(float));
// ar.readRawData((char *)m_vProperty,m_size*sizeof(float));
// END_READ_BLOCK( 1 );
// END_READ_OBJECT( ar );
// }
// return true;
// }
// bool DOUBLEPROPERTY::Serialize( CObjectArchive &ar )
// {
// if(ar.IsStoring())
// {
// BEGIN_WRITE_OBJECT( ar , 1 );
// BEGIN_WRITE_BLOCK( ar , 1 );
// ar << m_size;
// ar.writeRawData((char *)m_vProperty,m_size*sizeof(double));
// END_WRITE_BLOCK( ar , 1);
// END_WRITE_OBJECT( ar );
// }
// else
// {
// BEGIN_READ_OBJECT( ar,1 );
// BEGIN_READ_BLOCK( 1 );
// ar >> m_size;
// if(m_vProperty) delete m_vProperty;
// m_vProperty =new double[m_size];//(float *)malloc(sizeof(float));
// ar.readRawData((char *)m_vProperty,m_size*sizeof(double));
// END_READ_BLOCK( 1 );
// END_READ_OBJECT( ar );
// }
// }
void CreateDir(char *buf)
{
char buf1[MAX_PATH+1];
QDir adir;
if(!adir.mkdir(buf)) {
int j=0;
int i=0;
CString cs=buf;
for(i=2;i<strlen(buf);i++) {
if(buf[i]=='\\'||buf[i]=='/') {
if(!j) {
j++;
continue;
}
strcpy(buf1,buf);
buf1[i]='\0';
if(!adir.mkdir(buf1)) {
/*
if(GetLastError()!=ERROR_ALREADY_EXISTS) {
CString cs;
cs.Format("新建目录%s错误",buf1);
AfxMessageBox(cs);
}
*/
}
}
else if(i!=0&&i==strlen(buf)-1) {
QDir adir;
adir.mkdir(buf);
}
}
}
}
// CBaseFunApp
void PutScanDepthMes(char *Message,float stdep,float endep,char *outfilename)
{
int len=strlen(Message);
Message[len]='\0';
len+=1;
char buf[128];
sprintf(buf,"%.6f\0",stdep);
strcpy(&Message[len],buf);//跳过\0
len+=strlen(buf)+1;//buf+1
sprintf(buf,"%.6f\0",endep);
strcpy(&Message[len],buf);
len+=strlen(buf)+1;
if(outfilename) {
sprintf(buf,"%.3f\0",endep);
strcpy(&Message[len],outfilename);
}
else strcpy(&Message[len],"");
}
void GetTranMes(char *p,float *x,float *y,float *outsdep,float *outedep)
{
sscanf(p,"%f",x);
int len=strlen(p)+1;
if(len<=1)*x=-99999.;
p+=len;
sscanf(p,"%f",y);
len=strlen(p)+1;
if(len<=1)*y=-99999.;
p+=len;
sscanf(p,"%f",outsdep);
len=strlen(p)+1;
p+=len;
sscanf(p,"%f",outedep);
/*
if(*outsdep>*outedep)
{
float temp=*outsdep;
*outsdep=*outedep;
*outedep=temp;
}
*/
}
float pe_axp_(int *cc)
{
union {
float x1;
int x2;
}zz;
int a,b,c,c1,c2,c3;
short a1=0,b1=0;
short i,j;
zz.x2=*cc;
if(zz.x2==0)return zz.x1;
a=0x00000080&zz.x2;
b1=zz.x2&0x0000007f;
c=0xffffff00&zz.x2;
c1=(c>>8)&0xff;
i=128;
for(j=0;j<3;++j){
if(c1>=i)break;
i=i/2;
};
c1=(c&0x00007f00)<<8;
c2=(c&0x00ff0000)>>8;
c3=(c&0xff000000)>>24;
c=((c1|c2|c3)&0X00ffffff)<<j;
if(b1>=0x41){
b=(b1-0x41)*4+0x7f+3-j;
}
else
{
b=0x7e-(0x40-b1)*4-j;
}
zz.x2=0x80000000&(a<<24);
zz.x2=zz.x2|((b<<23)&0x7f800000);
zz.x2=zz.x2|(c&0x007fffff);
return zz.x1;
}
int int32_pe_axp_(int *number)
{
int num,hi1,hi2,lo1,lo2;
num=*number;
hi1=(num&0x000000ff)<<24;
hi2=(num&0x0000ff00)<<8;
lo1=(num&0x00ff0000)>>8;
lo2=(num&0xff000000)>>24;
return hi1|hi2|lo1|lo2;
}
short int int16_pe_axp_(short int *number)
{
short num,hi,lo;
num=*number;
hi=(num&0xff00)>>8;
lo=(num&0x00ff)<<8;
return hi|lo;
}
CString GetBinDir()
{
char bindir[MAX_PATH+1];
bindir[0]=0;
CString cs=GetBinDir(bindir);
return cs;
}
char *GetBinDir(char *str)
{
QString strPathTmp = QCoreApplication::applicationDirPath()+ QDir::separator() ;
QString strConfPath = QDir::toNativeSeparators( strPathTmp );
strcpy(str,strConfPath.toStdString().c_str());
return str;
}
bool InRange(float val,float min,float max)
{
if(min==-9999)
{
if(val<=max)return 1;
else return 0;
}
else if(max==-9999)
{
if(val>=min)return 1;
else return 0;
}
else
{
if(val>=min&&val<=max)return 1;
else return 0;
}
}
int GetPos(DATAFORMATMESSAGE ms,char *buf)
{
CString Type[20]={"C","B","I2S","I2","I4S","I4","FPS","FP","UI4"};
int DataTypeId,starPos=0,pos=0;
int NumCTerm[10];//各条件中枚举字符串的个数,如果不是"C"型NumCTerm[i]=0
CString Cstr[50][10];
for(int i=0;i<ms.TermNum;i++)
{
NumCTerm[i]=0;
if(strcmp(ms.m_Term[i].ValueType,"C")!=0)continue;
Cstr[i][0].Format("%s",ms.m_Term[i].ValueMin);
int n=1;
CString cmax;
cmax.Format("%s",ms.m_Term[i].ValueMax);
while(1)
{
char strTmp[5] = ".OR.";
int p=cmax.Find(strTmp);
if(p<=0)
{
//Cstr[i][n].Format("%s",cmax);
Cstr[i][n]=cmax;
n++;
break;
}
Cstr[i][n]=cmax.Left(p);
cmax=cmax.Right(cmax.GetLength()-p-4);
n++;
}
NumCTerm[i]=n;
}
while(1)
{
for(int i=0;i<ms.TermNum;i++)
{
char mbuf[5000];
int ist,iend;
if(ms.m_Term[i].SearchMode==0){ist=ms.m_Term[i].start-1;iend=ms.m_Term[i].end;}//strncpy(mbuf,&buf[ms.m_Term[i].start-1],ms.m_Term[i].end-ms.m_Term[i].start+1);//全范围搜索
else if(ms.m_Term[i].SearchMode==1){ist=ms.m_Term[i].start-1;iend=ms.m_Term[i].end;}//strncpy(mbuf,&buf[ms.m_Term[i].start-1],ms.m_Term[i].end-ms.m_Term[i].start+1);//定位搜索
else if(ms.m_Term[i].SearchMode==2){ist=starPos+ms.m_Term[i].start;iend=ms.m_Term[i].end;}//strncpy(mbuf,&buf[starPos+ms.m_Term[i].start-1],ms.m_Term[i].end);//在第一搜索基础上加偏移量搜索
for(int j=ist;j<iend;j++)mbuf[j-ist]=buf[j];
CString DataType;
DataType.Format("%s",ms.m_Term[i].ValueType);
DataTypeId=-1;
for(int j=0;j<20;j++)
{
if(DataType==Type[j])
{
DataTypeId=j;
break;
}
}
if(DataTypeId<0)return 0;
switch(DataTypeId)
{
case 0://C
{
for(int j=0;j<iend-ist;j++)
{
if(mbuf[j]==0)mbuf[j]=' ';
else if(mbuf[j]>='a'&&mbuf[j]<='z')
{
mbuf[j]-=32;
}
}
CString str1;
str1.Format("%s",mbuf);
int l=str1.GetLength();
for(int m=0;m<NumCTerm[i];m++)
{
pos=str1.Find(Cstr[i][m]);
if(pos>=0)break;
}
//pos=str1.Find(ms.m_Term[i].ValueMin);
//if(pos<0)pos=str1.Find(ms.m_Term[i].ValueMax);
if(pos<0)return 0;
break;
}
case 1://B
{
BYTE val=mbuf[0];
float min=atof(ms.m_Term[i].ValueMin);
float max=atof(ms.m_Term[i].ValueMax);
if(!InRange(val,min,max))goto loop1;
//if(val<min||val>max)goto loop1;
break;
}
case 2://I2S
{
short int val;
memcpy(&val,mbuf,2);
val=int16_pe_axp_(&val);
float min=atof(ms.m_Term[i].ValueMin);
float max=atof(ms.m_Term[i].ValueMax);
//if(val<min||val>max)goto loop1;
if(!InRange(val,min,max))goto loop1;
break;
}
case 3://I2
{
short int val;
memcpy(&val,mbuf,2);
int min=atoi(ms.m_Term[i].ValueMin);
int max=atoi(ms.m_Term[i].ValueMax);
//if(val<min||val>max)goto loop1;
if(!InRange(val,min,max))goto loop1;
break;
}
case 4://I4S
{
int val;
memcpy(&val,mbuf,2);
val=int32_pe_axp_(&val);
int min=atoi(ms.m_Term[i].ValueMin);
int max=atoi(ms.m_Term[i].ValueMax);
//if(val<min||val>max)goto loop1;
if(!InRange(val,min,max))goto loop1;
break;
}
case 5://I4
{
int val;
memcpy(&val,mbuf,2);
int min=atoi(ms.m_Term[i].ValueMin);
int max=atoi(ms.m_Term[i].ValueMax);
//if(val<min||val>max)goto loop1;
if(!InRange(val,min,max))goto loop1;
break;
}
case 6://FP
{
int ival;
float val;
memcpy(&ival,mbuf,4);
val=pe_axp_(&ival);
float min=atof(ms.m_Term[i].ValueMin);
float max=atof(ms.m_Term[i].ValueMax);
//if(val<min||val>max)goto loop1;
if(!InRange(val,min,max))goto loop1;
break;
}
case 7://FPPC
{
float val;
memcpy(&val,mbuf,4);
float min=atof(ms.m_Term[i].ValueMin);
float max=atof(ms.m_Term[i].ValueMax);
// if(val<min||val>max)goto loop1;
if(!InRange(val,min,max))goto loop1;
break;
}
case 8://unsigned long
{
unsigned long val;
memcpy(&val,mbuf,4);
float min=atof(ms.m_Term[i].ValueMin);
float max=atof(ms.m_Term[i].ValueMax);
//if(val<min||val>max)goto loop1;
if(!InRange(val,min,max))goto loop1;
break;
}
}
if(i==0)starPos=pos+ms.m_Term[i].start-1;
}
return starPos;
loop1: if(ms.m_Term[0].SearchMode==0)
{
CString ss;
ss.Format("%s",ms.m_Term[0].ValueMin);
int h=ss.GetLength();
ms.m_Term[0].start=starPos+1+h;
if(ms.m_Term[0].start>=ms.m_Term[0].end)return -1;
}
else return -1;
}
return starPos;
}
int GetBeginPosit(char *FileName,char* FormatName)
{
//读数据格式识别信息文件
char *pPath=new char[MAX_PATH+1];
GetBinDir(pPath);
CString csDir = pPath;
csDir+=CString("\\convertor\\");
CString csFile = csDir +"DataFormatMes.txt",str;
delete pPath;
FILE *fp=fopen(csFile.GetString(),"rt");
if(fp==NULL)
{
AfxMessageBox("打开数据格式识别信息文件DataFormatMes.txt出错");
return -1;
}
char strName[256],temp[256],buf[7][30];
fgets(strName,256,fp);
DATAFORMATMESSAGE ms;
int LenMax=0;
while ( !feof(fp) )
{
fgets(strName,256,fp);
CString str;
str.Format("%s",strName);
str.MakeUpper();
if(str.Find(FormatName,0)>=0)
{
sscanf(strName,"%s %d",temp,&ms.TermNum);
ms.FormatName.Format("%s",temp);
ms.FormatName.Delete(0,1);
ms.FormatName.Delete(ms.FormatName.GetLength()-1,1);
ms.m_Term=new TERM[ms.TermNum];
for(int i=0;i<ms.TermNum;i++)
{
fgets(strName,256,fp);
sscanf(strName,"%s%s%s%s%s%s%s",ms.m_Term[i].characterization,buf[0],buf[1],buf[2],
ms.m_Term[i].ValueType,ms.m_Term[i].ValueMin,ms.m_Term[i].ValueMax);//buf[3],buf[4],buf[5],buf[6]);//%d %d %d %d %s %s %s",temp,&ms.m_Term[i].SearchMode,&ms.m_Term[i].start,
// &ms.m_Term[i].end,&buf[0],&buf[1],&buf[2]);
sscanf(buf[0],"%d",&ms.m_Term[i].SearchMode);
sscanf(buf[1],"%d",&ms.m_Term[i].start);
sscanf(buf[2],"%d",&ms.m_Term[i].end);
if(LenMax<ms.m_Term[i].end)LenMax=ms.m_Term[i].end;
}
break;
}
}
fclose(fp);
//格式识别
char *buffer;
buffer=new char[LenMax+1];
fp=fopen(FileName,"rb");
fread(buffer,1,LenMax,fp);
fclose(fp);
buffer[LenMax]=0;
int pos=GetPos(ms,buffer);
delete ms.m_Term;
delete buffer;
return pos;
}
CString GetRegDir(char *Name)
{
// TODO: Add your control notification handler code here
char Path[MAX_PATH+1];
CString strDLL=GetBinDir(Path);
if(strDLL.ReverseFind('\\')>-1) {
strDLL=strDLL.Left(strDLL.ReverseFind('\\'));
strcpy(Path,strDLL.GetString());
}
strcat(Path,"\\");
strcat(Path,Name);
return Path;
}
bool HaveSameCurve(int NumLog,int *OutCurveNo,char **OutCurve)
{
if(NumLog>1024) return 1;
for(int i=0;i<NumLog-1;i++)
{
char temp[1000];
strcpy(temp,OutCurve[i]);
if(strlen(temp)==0&&OutCurveNo[i]>=0) {
AfxMessageBox("输出曲线名不能为空!");
return 1;
break;
}
for(int j=i+1;j<NumLog;j++)
{
if((strcmp(OutCurve[i],OutCurve[j])==0)&&(OutCurveNo[i]>=0))
{
strcat(temp,"_1");
delete OutCurve[i];
OutCurve[i]=new char[strlen(temp)+1];
strcpy(OutCurve[i],temp);
}
}
}
/*
CString str,mes="";
int SameNameCurveNum=0;
for(int i=0;i<NumLog-1;i++)
{
//SameNameCurveNum=0;
for(int j=i+1;j<NumLog;j++)
{
if((strcmp(OutCurve[i],OutCurve[j])==0)&&(OutCurveNo[i]>=0))
{
//OutCurveNo[j]=-1;
sprintf(OutCurve[j],"%s%d",OutCurve[j],SameNameCurveNum);
str.Format("第%d条曲线和第%d条曲线的输出曲线名相同,输出名字为%s\r\n",i,j,OutCurve[i]);
mes+=str;
SameNameCurveNum++;
break;
}
}
}
if(SameNameCurveNum>0)
{
for(int i=0;i<NumLog-1;i++)
{
char temp[64];
strcpy(temp,OutCurve[i]);
for(int j=i+1;j<NumLog;j++)
{
if((strcmp(OutCurve[i],OutCurve[j])==0)&&(OutCurveNo[i]>=0))
{
strcat(temp,"+");
strcpy(OutCurve[j],temp);
}
}
}
}*/
return 0;
}
void fft842(int INN, int N, float *X, float *Y)
{
int I,IJ,JI,IPASS,IW,J,J1,J2,J3,J4,J5,J6,J7,J8,J9,J10,J11,J12,J13,J14,
LENGT,M,M1,M2,M3,M4,N2POW,N8POW,NT,NTHPO,NXTLT,L[15];
float FI,FN,R;
float XX[1025],YY[1025];
IW=0;
// PI2=8.*atan(1.);
// PI7=1./sqrt(2.);
for(I=1;I<=15;I++)
{
M=I;
NT=pow(2.0,I);
if(N == NT) goto L_20;
}
// exit(0);
L_20:
N2POW=M;
NTHPO=N;
FN=1.0f*NTHPO;
if(INN != 1)
{
for(I=0;I<NTHPO;I++)
Y[I]=-Y[I];
}
N8POW=N2POW/3;
if(N8POW != 0)
{
for(IPASS=0;IPASS<N8POW;IPASS++)
{
NXTLT=pow(2.0,N2POW-3*(IPASS+1));
LENGT=8*NXTLT;
r8tx(NXTLT,NTHPO,LENGT,X,Y);//
//X+NXTLT,X+NXTLT*2,X+NXTLT*3,X+NXTLT*4,X+NXTLT*5,X+NXTLT*6,
// X+NXTLT*7,Y,Y+NXTLT,Y+NXTLT*2,Y+NXTLT*3,Y+NXTLT*4,Y+NXTLT*5,Y+NXTLT*6,Y+NXTLT*7);
}
}
M1=N2POW-3*N8POW-1;
if(M1==0)
r2tx(NTHPO,X,Y);//X+1,Y,Y+1);
else if(M1>0)
r4tx(NTHPO,X,Y);//X+1,X+2,X+3,Y,Y+1,Y+2,Y+3);
for(J=1;J<=15;J++)
{
L[J-1]=1;
M2=J-N2POW;
if(M2 <=0 )
{
M3=N2POW-J+1;
L[J-1]=pow(2.0,M3);
}
}
IJ=1;
for(J1=1;J1<=N;J1++)
{
XX[J1]=X[J1-1];
YY[J1]=Y[J1-1];
}
for(J1=1;J1<=L[14];J1++)
{
for(J2=J1;J2<=L[13];J2+=L[14])
{
for(J3=J2;J3<=L[12];J3+=L[13])
{
for(J4=J3;J4<=L[11];J4+=L[12])
{
for(J5=J4;J5<=L[10];J5+=L[11])
{
for(J6=J5;J6<=L[9];J6+=L[10])
{
for(J7=J6;J7<=L[8];J7+=L[9])
{
for(J8=J7;J8<=L[7];J8+=L[8])
{
for(J9=J8;J9<=L[6];J9+=L[7])
{
for(J10=J9;J10<=L[5];J10+=L[6])
{
for(J11=J10;J11<=L[4];J11+=L[5])
{
for(J12=J11;J12<=L[3];J12+=L[4])
{
for(J13=J12;J13<=L[2];J13+=L[3])
{
for(J14=J13;J14<=L[1];J14+=L[2])
{
for(JI=J14;JI<=L[0];JI+=L[1])
{
M4=IJ-JI;
if(M4<0)
{
R=XX[IJ];
XX[IJ]=XX[JI];
XX[JI]=R;
FI=YY[IJ];
YY[IJ]=YY[JI];
YY[JI]=FI;
}
IJ++;
}
}
}
}
}
}
}
}
}
}
}
}
}
}
}
//
for(J1=1;J1<=N;J1++)
{
X[J1-1]=XX[J1];
Y[J1-1]=YY[J1];
}
if(INN == 1)
{
for(I=0;I<NTHPO;I++)
{
X[I]=X[I]/FN;
Y[I]=Y[I]/FN;
}
}
else
{
for(I=0;I<NTHPO;I++)
Y[I]=-Y[I];
}
}
void r2tx(int NTHPO, float *CR0, float *CI0)//1, float *CI0, float *CI1)
{
int K;
float FI1,R1;
for(K=0;K<NTHPO;K+=2)
{
R1=CR0[K]+CR0[K+1];
CR0[K+1]=CR0[K]-CR0[K+1];
CR0[K]=R1;
FI1=CI0[K]+CI0[K+1];
CI0[K+1]=CI0[K]-CI0[K+1];
CI0[K]=FI1;
}
return;
}
void r4tx(int NTHPO, float *CR0, float *CI0)//, float *CR2, float *CR3, float *CI0, float *CI1, float *CI2, float *CI3)
{
int K;
float FI1,FI2,FI3,FI4,R1,R2,R3,R4;
for(K=0;K<NTHPO;K+=4)
{
R1=CR0[K]+CR0[K+2];
R2=CR0[K]-CR0[K+2];
R3=CR0[K+1]+CR0[K+3];
R4=CR0[K+1]-CR0[K+3];
FI1=CI0[K]+CI0[K+2];
FI2=CI0[K]-CI0[K+2];
FI3=CI0[K+1]+CI0[K+3];
FI4=CI0[K+1]-CI0[K+3];
CR0[K]=R1+R3;
CI0[K]=FI1+FI3;
CR0[K+1]=R1-R3;
CI0[K+1]=FI1-FI3;
CR0[K+2]=R2-FI4;
CI0[K+2]=FI2+R4;
CR0[K+3]=R2+FI4;
CI0[K+3]=FI2-R4;
}
}
void r8tx(int NXTLT,int NTHPO,int LENGT,float *CR0,float *CI0)
{
float PI2=8.*atan(1.);
float PI7=1./sqrt(2.);
int J,K;
float AI0,AI1,AI2,AI3,AI4,AI5,AI6,AI7,
AR0,AR1,AR2,AR3,AR4,AR5,AR6,AR7,
BI0,BI1,BI2,BI3,BI4,BI5,BI6,BI7,
BR0,BR1,BR2,BR3,BR4, BR5,BR6,BR7,
C1,C2,C3,C4,C5,C6,C7,S1,S2,S3,S4,S5,S6,S7,ARG,SCALE,TI,TR;
SCALE=PI2/LENGT;
for(J=0;J<NXTLT;J++)
{
ARG=J*SCALE;
C1=cos(ARG);
S1=sin(ARG);
C2=C1*C1-S1*S1;
S2=C1*S1+C1*S1;
C3=C1*C2-S1*S2;
S3=C2*S1+S2*C1;
C4=C2*C2-S2*S2;
S4=C2*S2+C2*S2;
C5=C2*C3-S2*S3;
S5=C3*S2+S3*C2;
C6=C3*C3-S3*S3;
S6=C3*S3+C3*S3;
C7=C3*C4-S3*S4;
S7=C4*S3+S4*C3;
for(K=J;K<NTHPO;K=K+LENGT)
{
AR0=CR0[K]+CR0[K+4*NXTLT];
AR1=CR0[K+NXTLT]+CR0[K+5*NXTLT];
AR2=CR0[K+2*NXTLT]+CR0[K+6*NXTLT];
AR3=CR0[K+3*NXTLT]+CR0[K+7*NXTLT];
AR4=CR0[K]-CR0[K+4*NXTLT];
AR5=CR0[K+NXTLT]-CR0[K+5*NXTLT];
AR6=CR0[K+2*NXTLT]-CR0[K+6*NXTLT];
AR7=CR0[K+3*NXTLT]-CR0[K+7*NXTLT];
AI0=CI0[K]+CI0[K+4*NXTLT];
AI1=CI0[K+NXTLT]+CI0[K+5*NXTLT];
AI2=CI0[K+2*NXTLT]+CI0[K+6*NXTLT];
AI3=CI0[K+3*NXTLT]+CI0[K+7*NXTLT];
AI4=CI0[K]-CI0[K+4*NXTLT];
AI5=CI0[K+NXTLT]-CI0[K+5*NXTLT];
AI6=CI0[K+2*NXTLT]-CI0[K+6*NXTLT];
AI7=CI0[K+3*NXTLT]-CI0[K+7*NXTLT];
BR0=AR0+AR2;
BR1=AR1+AR3;
BR2=AR0-AR2;
BR3=AR1-AR3;
BR4=AR4-AI6;
BR5=AR5-AI7;
BR6=AR4+AI6;
BR7=AR5+AI7;
BI0=AI0+AI2;
BI1=AI1+AI3;
BI2=AI0-AI2;
BI3=AI1-AI3;
BI4=AI4+AR6;
BI5=AI5+AR7;
BI6=AI4-AR6;
BI7=AI5-AR7;
CR0[K]=BR0+BR1;
CI0[K]=BI0+BI1;
if(J > 0)
{
CR0[K+NXTLT]=C4*(BR0-BR1)-S4*(BI0-BI1);
CI0[K+NXTLT]=C4*(BI0-BI1)+S4*(BR0-BR1);
CR0[K+2*NXTLT]=C2*(BR2-BI3)-S2*(BI2+BR3);
CI0[K+2*NXTLT]=C2*(BI2+BR3)+S2*(BR2-BI3);
CR0[K+3*NXTLT]=C6*(BR2+BI3)-S6*(BI2-BR3);
CI0[K+3*NXTLT]=C6*(BI2-BR3)+S6*(BR2+BI3);
TR=PI7*(BR5-BI5);
TI=PI7*(BR5+BI5);
CR0[K+4*NXTLT]=C1*(BR4+TR)-S1*(BI4+TI);
CI0[K+4*NXTLT]=C1*(BI4+TI)+S1*(BR4+TR);
CR0[K+5*NXTLT]=C5*(BR4-TR)-S5*(BI4-TI);
CI0[K+5*NXTLT]=C5*(BI4-TI)+S5*(BR4-TR);
TR=-PI7*(BR7+BI7);
TI=PI7*(BR7-BI7);
CR0[K+6*NXTLT]=C3*(BR6+TR)-S3*(BI6+TI);
CI0[K+6*NXTLT]=C3*(BI6+TI)+S3*(BR6+TR);
CR0[K+7*NXTLT]=C7*(BR6-TR)-S7*(BI6-TI);
CI0[K+7*NXTLT]=C7*(BI6-TI)+S7*(BR6-TR);
}
else
{
CR0[K+NXTLT]=BR0-BR1;
CI0[K+NXTLT]=BI0-BI1;
CR0[K+2*NXTLT]=BR2-BI3;
CI0[K+2*NXTLT]=BI2+BR3;
CR0[K+3*NXTLT]=BR2+BI3;
CI0[K+3*NXTLT]=BI2-BR3;
TR=PI7*(BR5-BI5);
TI=PI7*(BR5+BI5);
CR0[K+4*NXTLT]=BR4+TR;
CI0[K+4*NXTLT]=BI4+TI;
CR0[K+5*NXTLT]=BR4-TR;
CI0[K+5*NXTLT]=BI4-TI;
TR=-PI7*(BR7+BI7);
TI=PI7*(BR7-BI7);
CR0[K+6*NXTLT]=BR6+TR;
CI0[K+6*NXTLT]=BI6+TI;
CR0[K+7*NXTLT]=BR6-TR;
CI0[K+7*NXTLT]=BI6-TI;
}
}
}
}
bool IsNum(CString str)
{
for(int i=0;i<str.GetLength();i++)
{
char c=str.GetAt(i);
if((c<='0'||c>='9')&&c!='.')
return false;
}
return true;
}
bool MiddleSmooth(float* _SourceData,int Count,int WindowLen)
{
using namespace std;
if (_SourceData == NULL) return false;
if (WindowLen < 2) return false;
if (Count < WindowLen) return false;
int DataCount = Count;
float* RetData = new float[DataCount];
for(int k=0;k<DataCount;k++)RetData[k]=_SourceData[k];
for (int i = WindowLen / 2; i < DataCount - WindowLen / 2; i++)
{
list <float> l;
for (int k = i - WindowLen / 2; k < i + WindowLen / 2 + 1; k++)l.push_back(RetData[k]);
l.sort();
for(int kk=0;kk<WindowLen / 2;kk++)l.pop_back();
_SourceData[i] =l.back();
}
delete RetData;
return true;
}
QStringList GetStringList(QByteArray line,bool IsSpa,bool IsTab,bool IsCom,bool IsSem,bool DelDubSpa)
{
QTextCodec *codec = QTextCodec::codecForName("UTF8");
if(!IsCom)line.replace(",","");//如果逗号不是分隔符,清除逗号
if(!IsSpa)line.replace(" ","");
else
{
//连续空格视为单个处理
if(DelDubSpa)while(line.indexOf(" ")>=0)line.replace(" "," ");
// line=line.simplified();
line.replace(" ",",");
}
if(!IsTab)line.replace(" ","");//tab
else line.replace(" ",",");
if(!IsSem)line.replace(";","");
else line.replace(";",",");
// line.remove("\r");
// line.remove("\n");
QList<QByteArray> list=line.split(',');
QStringList List;
for(int i=0;i<list.size();i++)
{
if(list[i].size()) {
unsigned char c=list[i].at(0)&0xe0;
if(c==0xe0)
{
List.append(codec->toUnicode(list[i],list[i].size()));
}
else List.append(list[i]);
}
}
return List;
}
void valueperlayers(float *csdata,int cs,float *grm,float*dep,int mpoint,float *resc,float*res,float gryz,int isDelerr)
{
float*grmlb=new float[mpoint];
int lbn,n,aa,bb;
int i,j,k;
//!异常值剔除
if(isDelerr) for(i=1;i<mpoint;i++)
{
if(grm[i]>1000.0||grm[i]<=20.0) grm[i]=grm[i-1];
}
//!光滑滤波(去除毛刺)
for(int i=0;i<mpoint;i++) grmlb[i]=grm[i];
lbn=10;
float yz=10.0;
for(i=lbn;i<mpoint-lbn;i++)
{
double sum=0;
for(int j=i-lbn-1;j<i+lbn;j++){
sum+=grm[j];
}
float temp=sum/(2.0*lbn+1);
if(grm[i]>=temp+yz||grm[i]<=temp-yz) grmlb[i]=temp;
}
lbn=2;
res=grmlb;
for(i=lbn;i<mpoint-lbn;i++)
{
double sum=0;
for(int j=i-lbn-1;j<i+lbn;j++){
sum+=grmlb[j];
}
res[i]=sum/(2.0*lbn+1.0);
}
for(j=0;j<cs;j++){
aa=0;
bb=0;
for(i=1;i<mpoint-1;i++){
if(csdata[j]>=dep[i]&&csdata[j]<=dep[i+1]) aa=i;
if(csdata[j+1]>=dep[i]&&csdata[j+1]<=dep[i+1]) bb=i;
}
float temp=0;
n=0;
for(k=aa;k<bb;k++){
n=n+1;
temp=temp+res[k];
}
temp=res[(int)(aa+bb)/2];// !取层中间点
resc[j]=temp;
for(k=aa-1;k<bb-1;k++){
res[k]=temp;
}
}
for(i=1;i<mpoint;i++){
if(fabs(res[i]-res[i-1])<=gryz) {
res[i]=res[i-1];
}
}
for(j=1;j<cs;j++){
if(fabs(resc[j]-resc[j-1])<=gryz){
resc[j]=resc[j-1];
}
}
return;
}
float* AutoCreateLayers(float *ldepth,float*grm,int lnsamp,float minceng,int &count,float gryz,int n,int isdellerr)
{
float *lwdgr=new float[lnsamp];
float *slwd=new float[lnsamp];
float *ratio=new float[lnsamp];
float *gedge=new float[lnsamp];
float *gratio=new float[lnsamp];
float *dd=new float[lnsamp];
float *gpratio=new float[lnsamp];
float *cengph=new float[lnsamp];
float *edge=new float[lnsamp];
float *cenggr=new float[lnsamp];
float *pratio=new float[lnsamp];
float *zceng=new float[lnsamp];
for(int i=0;i<lnsamp;i++)
{
slwd[i]=grm[i];
}
float ldet=ldepth[2]-ldepth[1];
// !按视电阻率分层
// int n=10;// !窗长
for(int i=n/2;i<lnsamp-n/2-1;i++)
{
gratio[i]=0.0;
for(int j=1-1;j<n/2;j++)
{
gratio[i]=gratio[i]+slwd[j+i]-slwd[i-j];
}
gratio[i]=fabs(gratio[i])/n/(n+1)/ldet;
dd[i]=ldepth[i];
}
int mm=0;
int bbrla=1;
gedge[0]=ldepth[0];
for(int i=n/2;i<lnsamp-n/2-1;i++)
{
mm=mm+1;
gedge[mm]=ldepth[i];
gpratio[mm]=gratio[i];
gpratio[mm+1]=gratio[i+1];
if((gedge[mm]-gedge[0])>=minceng)
{
if(mm>=1){
if(((gpratio[mm-1]+0.0001)<=gpratio[mm])&&(gpratio[mm]>=(gpratio[mm+1]+0.0001)))
{
bbrla=bbrla+1;
cengph[bbrla-1]=gedge[mm];
cengph[bbrla]=gedge[mm]+0.01;
gedge[0]=gedge[mm];
}
}
}
}
cengph[0]=ldepth[0];
cengph[bbrla]=ldepth[lnsamp-1];
/*
//!按视自然伽马曲线gr分层
n=3;
for(int i=n/2;i<lnsamp-n/2-1;i++)
{
ratio[i]=0.0;
for(int j=0;j<n/2;j++)
{
ratio[i]=ratio[i]+lwdgr[j+i]-lwdgr[i-j];
}
ratio[i]=fabs(ratio[i])/n/(n+1)/ldet;
dd[i]=ldepth[i];
}
mm=0;
int bbgr=1;
edge[0]=ldepth[0];
for(int i=n/2;i<lnsamp-n/2-1;i++){
mm=mm+1;
edge[mm]=ldepth[i];
pratio[mm]=ratio[i];
pratio[mm+1]=ratio[i+1];
if((edge[mm]-edge[0])>=minceng){
if(mm>=1){
if(((pratio[mm-1]+0.0001)<=pratio[mm])&&(pratio[mm]>=(pratio[mm+1]+0.0001))){
bbgr=bbgr+1;
cenggr[bbgr-1]=edge[mm];
cenggr[bbgr]=edge[mm]+0.01;
edge[0]=edge[mm];
}
}
}
}
cenggr[0]=ldepth[0];
cenggr[bbgr]=ldepth[lnsamp-1];
count=-1;
zceng[0]=cengph[0];
for(int i=0;i<bbrla-1;i++)
{
count=count+1;
zceng[count]=cengph[i];
for(int ii=0;ii<bbgr-1;ii++)
{
if(((cenggr[ii]-cengph[i])>=minceng)
&&((cengph[i+1]-cenggr[ii])>=minceng)
&&((cenggr[ii]-cenggr[ii-1])>=minceng))
{
count=count+1;
zceng[count]=cenggr[ii];
}
}
}
zceng[count+1]=cengph[bbrla];
*/
count=-1;
zceng[0]=cengph[0];
for(int i=0;i<bbrla-1;i++)
{
count=count+1;
zceng[i]=cengph[i];
}
// count=count+1;
if(cengph[bbrla-1]>zceng[bbrla-2])zceng[count+1]=cengph[bbrla-1];
delete slwd;
delete gedge;
delete gratio;
delete dd;
delete gpratio;
delete cengph;
delete edge;
delete cenggr;
delete pratio;
if(gryz == 0) return zceng;
float* resc = new float[count+1];
float* res = new float[lnsamp];
if(gryz == 0)
gryz = 5;
float *nceng=new float[count+1];
valueperlayers(zceng,count-1,grm,ldepth,lnsamp,resc,res,gryz,isdellerr);
int ncount=0;
for(int i =0;i<count+1;i++)
{
if(resc[i]==resc[i+1]) continue;
nceng[ncount++] = zceng[i];
}
if(nceng[ncount-1]!=zceng[count+1])
{
nceng[ncount++] = zceng[count+1];
}
delete[] resc;
delete[] res;
delete zceng;
zceng=nceng;
count=ncount;
return zceng;
}
int CalculateDataNo(double* data1,int nsamp1,double data)
{
if(data >= data1[nsamp1-1])
return nsamp1-1;
if(data <= data1[0])
return 0;
int No = -1;
double min = 9999;
for(int j = 0;j<nsamp1;j++)
{
if(data1[j] - data > 0 && data1[j] - data < min)
{
min = data1[j] - data;
No = j;
//break;
}
}
return No;
}
int CalculateDataNo(float* data1,int nsamp1,double data)
{
if(data >= data1[nsamp1-1])
return nsamp1-1;
if(data <= data1[0])
return 0;
int No = -1;
double min = 9999;
for(int j = 0;j<nsamp1;j++)
{
if(data1[j] - data > 0 && data1[j] - data < min)
{
min = data1[j] - data;
No = j;
//break;
}
}
return No;
}
void CalculateLayer(double* DepthData,double* ResData,int nsamp,char* LayerPath,int LayerPathLeng,float MinVAL,float MaxVAL,float minLayerThick)
{
typedef void(*LAYERSCOMPUTE)(double*,double*,int,double,char*,int);
char bin[256];
QString bin1 = GetBinDir(bin);
bin1 += "\\app\\GetEveryLayerFillValue.dll";
QLibrary lib(bin1);
if (!lib.load())
{
QMessageBox::information(NULL,"","加载GetEveryLayerFillValue.dll失败");
return;
}
LAYERSCOMPUTE lc=(LAYERSCOMPUTE)lib.resolve("LAYERSCOMPUTE");
if (!lc)
{
QMessageBox::information(NULL,"","获取LAYERSCOMPUTE函数失败");
return;
}
if(MinVAL > MaxVAL)
{
float temp = MaxVAL;
MaxVAL = MinVAL;
MinVAL = temp;
}
int staNo,endNo;
staNo = CalculateDataNo(DepthData,nsamp,MinVAL);
endNo = CalculateDataNo(DepthData,nsamp,MaxVAL);
nsamp = endNo-staNo;
double* DepthData1 = new double[nsamp];
double* ResData1 = new double[nsamp];
for(int i = staNo;i<staNo+nsamp;i++)
{
DepthData1[i-staNo] = DepthData[i];
ResData1[i-staNo] = ResData[i];
}
double minceng =minLayerThick;//最小层厚
if(minceng == 0)
minceng = 0.5;
lc(DepthData1,ResData1,nsamp,minceng,LayerPath,LayerPathLeng);
delete[] DepthData1;
delete[] ResData1;
}
void CalculateLayer(float* DepthData,float* ResData,int nsamp,char* LayerPath,int LayerPathLeng,float MinVAL,float MaxVAL,float minLayerThick)
{
typedef void(*LAYERSCOMPUTE)(double*,double*,int,double,char*,int);
char bin[256];
QString bin1 = GetBinDir(bin);
bin1 += "\\app\\GetEveryLayerFillValue.dll";
QLibrary lib(bin1);
if (!lib.load())
{
QMessageBox::information(NULL,"","加载GetEveryLayerFillValue.dll失败");
return;
}
LAYERSCOMPUTE lc=(LAYERSCOMPUTE)lib.resolve("LAYERSCOMPUTE");
if (!lc)
{
QMessageBox::information(NULL,"","获取LAYERSCOMPUTE函数失败");
return;
}
if(MinVAL > MaxVAL)
{
float temp = MaxVAL;
MaxVAL = MinVAL;
MinVAL = temp;
}
int staNo,endNo;
staNo = CalculateDataNo(DepthData,nsamp,MinVAL);
endNo = CalculateDataNo(DepthData,nsamp,MaxVAL);
nsamp = endNo-staNo;
if(nsamp<1) return;
double* DepthData1 = new double[nsamp];
double* ResData1 = new double[nsamp];
for(int i = staNo;i<staNo+nsamp;i++)
{
DepthData1[i-staNo] = DepthData[i];
ResData1[i-staNo] = ResData[i];
}
double minceng =minLayerThick;//最小层厚
if(minceng == 0)
minceng = 0.5;
lc(DepthData1,ResData1,nsamp,minceng,LayerPath,LayerPathLeng);
delete[] DepthData1;
delete[] ResData1;
}
void InstValue(float outsdep,float enddepth,std::vector<float> &md,std::vector<float> &_AryCurve,std::vector<float>&curAry,float rlev,int DataPoint,bool isAngle,bool isJG)
{
int size = (enddepth - outsdep) / rlev + 1.5 + 1;
curAry.clear();
float curdepth=outsdep;
curAry.reserve(size*DataPoint);
int index=0;
int indf=0;
int inde=0;
float depthf=curdepth;
float depthe=enddepth;
DWORD point=0;
while (curdepth<=enddepth)
{
if(index==md.size()-1)
{
if(curdepth==enddepth)
{
for(int n=0;n<DataPoint;n++)
{
float v=_AryCurve[index*DataPoint+n];
if(isJG) {
v=100;
}
curAry.push_back(v);
}
point++;
curdepth=outsdep+point*rlev;
continue;
}
while(curdepth<=enddepth)
{
float v=-9999.0;
if(isJG) {
v=0;
}
curAry.push_back(v);
point++;
curdepth=outsdep+point*rlev;
}
break;
}
float depth1=md[index];
if(fabsf(curdepth-depth1)>=rlev) {
index++;
if (index>=md.size())
break;
if(md[index-1] < curdepth && md[index] > curdepth)//后步过界,回退
index--;
else if(curdepth<md[index-1]){
}
else {
if (index>=md.size())
break;
continue;
}
}
if(index+1>=md.size()) break;
float depth2=md[index+1];
float rate=0;
if(depth2==depth1) {
for(int n=0;n<DataPoint;n++)
{
float val0=_AryCurve[index*DataPoint+n];
curAry.push_back(val0);
}
}
else {
rate=(curdepth-depth1)/(depth2-depth1);
for(int n=0;n<DataPoint;n++)
{
float val0=_AryCurve[index*DataPoint+n];
float val1=_AryCurve[(index+1)*DataPoint+n];
if(val0==-9999||val1==-9999||val1==-999.25||val0==-999.25)
{
float v=-9999.0;
if(isJG) {
if(curdepth==depth1||curdepth==depth2) {
if(curdepth<depth1) v=-100;
else if(curdepth>=depth1) v=100;
}
else if(fabs((int)(curdepth*1000+0.5)/1000.0-(int)(depth1*1000+0.5)/1000.0)<rlev-0.0001)
{
if(curdepth<depth1) v=-100;
else if(curdepth>=depth1) v=100;
}
else if(fabs((int)(curdepth*1000+0.5)/1000.0-(int)(depth2*1000+0.5)/1000.0)<rlev-0.0001)
{
if(curdepth<depth1) v=-100;
else if(curdepth>=depth1) v=100;
}
else v=0;
}
curAry.push_back(v);
}
else
{
if(isAngle) {
if(curdepth<depth2) rate=0;
else rate=1;
}
float v=val0+rate*(val1-val0);
if(isJG) {
if(fabs(curdepth-depth1)<rlev)
{
if(curdepth<depth1) v=-100;
else if(curdepth>=depth1) v=100;
}
else v=0;
}
curAry.push_back(v);
}
}
}
if(curdepth > enddepth && curdepth < enddepth+rlev){
for(int n=0;n<DataPoint;n++)
{
float v = _AryCurve[index*DataPoint+n];
if(isJG) curAry.push_back(0);
else curAry.push_back(v);
}
}
point++;
curdepth=outsdep+point*rlev;
while(curdepth>=depth2) {
index++;
if(index+1>=md.size()) break;
depth2=md[index+1];
}
if(index>=md.size()) break;
}
int count=curAry.size()/DataPoint;
if(count<size)
{
for(int j=count;j<size;j++)
{
for(int n=0;n<DataPoint;n++)
{
float v=-9999.0;
if(isJG) {
v=0;
}
curAry.push_back(v);
}
}
}
}
float Filt(float* b, float* Fil)
{
float num = 0.0f;
int index = 0;
do
{
if (index == 0)
num += b[5] * Fil[index];
else
num = (float) ((double) num + (double) b[ (5 + index)] * (double) Fil[index] + (double) b[ (5 - index)] * (double) Fil[index]);
{ ++index; }
}
while (index <= 5);
return num;
}
void FilterCurve(
float* inpcurve,
enumFilterType Ftype,
int numrec,
float* fUDef,
float Empty)
{
float* numArray1 = new float[6];
float* numArray2 = new float[6];
float* numArray3 = new float[6];
float* numArray4 = new float[6];
float* numArray5 = new float[6];
float* b = new float[11];
float* numArray6 = new float[6];
numArray1[0] = 0.3292f;
numArray1[1] = 0.254f;
numArray1[2] = 0.1021f;
numArray1[3] = -0.0073f;
numArray1[4] = -0.017f;
numArray1[5] = 0.0036f;
numArray2[0] = 3.0f / 16.0f;
numArray2[1] = 0.1704f;
numArray2[2] = 0.1266f;
numArray2[3] = 0.0735f;
numArray2[4] = 0.0298f;
numArray2[5] = 0.0059f;
numArray3[0] = 0.14272f;
numArray3[1] = 0.1217f;
numArray3[2] = 0.11401f;
numArray3[3] = 0.09303f;
numArray3[4] = 0.06583f;
numArray3[5] = 89.0f / (961.0f * 2.71828);
float* numArray7 = fUDef;
numArray5[0] = 1.0f;
numArray5[1] = 0.0f;
numArray5[2] = 0.0f;
numArray5[3] = 0.0f;
numArray5[4] = 0.0f;
numArray5[5] = 0.0f;
float* numArray8 = new float[ (numrec + 1)];
float* Fil;
switch (Ftype)
{
case enumFilterType::fltWeak:
Fil = numArray1;
break;
case enumFilterType::fltNormal:
Fil = numArray2;
break;
case enumFilterType::fltStrong:
Fil = numArray3;
break;
case enumFilterType::fltUDef:
Fil = numArray7;
break;
case enumFilterType::fltNone:
Fil = numArray5;
break;
default:
int num1 = 0;//(int) Interaction.MsgBox((object) "Wrong Filter Type! Filter Set to Weak!");
Fil = numArray1;
break;
}
int num2 = numrec;
int index1 = 1;
while (index1 <= num2)
{
if ((double) inpcurve[index1] == (double) Empty)
{
numArray8[index1] = Empty;
}
else
{
int index2 = -1;
int num3 = (index1 - 5);
int num4 = (index1 + 5);
int index3 = num3;
while (index3 <= num4)
{
{ ++index2; }
b[index2] = !(index3 < 1 | index3 > numrec) ? ((double) inpcurve[index3] != (double) Empty ? inpcurve[index3] : inpcurve[index1]) : inpcurve[index1];
{ ++index3; }
}
numArray8[index1] = Filt(b, Fil);
}
{ ++index1; }
}
for(int i=0;i<numrec;i++) inpcurve[i] = numArray8[i+1];
delete numArray1;
delete numArray2;
delete numArray3;
delete numArray4;
delete numArray5;
delete b;
delete numArray6;
delete numArray8;
}
float FilterValue(
float* inpcurve,
enumFilterType Ftype,
float* fUDef,
float Empty)
{
float* numArray1 = new float[6];
float* numArray2 = new float[6];
float* numArray3 = new float[6];
float* numArray4 = new float[6];
float* numArray5 = new float[6];
float* numArray6 = new float[11];
float* numArray7 = new float[6];
numArray1[0] = 0.3292f;
numArray1[1] = 0.254f;
numArray1[2] = 0.1021f;
numArray1[3] = -0.0073f;
numArray1[4] = -0.017f;
numArray1[5] = 0.0036f;
numArray2[0] = 3.0f / 16.0f;
numArray2[1] = 0.1704f;
numArray2[2] = 0.1266f;
numArray2[3] = 0.0735f;
numArray2[4] = 0.0298f;
numArray2[5] = 0.0059f;
numArray3[0] = 0.14272f;
numArray3[1] = 0.1217f;
numArray3[2] = 0.11401f;
numArray3[3] = 0.09303f;
numArray3[4] = 0.06583f;
numArray3[5] = 89.0f / (961.0f * 2.71828);
float* numArray8 = fUDef;
numArray5[0] = 1.0f;
numArray5[1] = 0.0f;
numArray5[2] = 0.0f;
numArray5[3] = 0.0f;
numArray5[4] = 0.0f;
numArray5[5] = 0.0f;
float* Fil;
switch (Ftype)
{
case enumFilterType::fltWeak:
Fil = numArray1;
break;
case enumFilterType::fltNormal:
Fil = numArray2;
break;
case enumFilterType::fltStrong:
Fil = numArray3;
break;
case enumFilterType::fltUDef:
Fil = numArray8;
break;
case enumFilterType::fltNone:
Fil = numArray5;
break;
default:
int num = 0;//(int) Interaction.MsgBox((object) "Wrong Filter Type! Filter Set to Weak!");
Fil = numArray1;
break;
}
return Filt(inpcurve, Fil);
delete numArray1;
delete numArray2;
delete numArray3;
delete numArray4;
delete numArray5;
delete numArray6;
delete numArray7;
}
void FilterMatrixV(
float** inpcurve,
enumFilterType Ftype,
int numrec,
int numChannel,
float* fUDef,
float Empty)
{
float* numArray1 = new float[6];
float* numArray2 = new float[6];
float* numArray3 = new float[6];
float* numArray4 = new float[6];
float* b = new float[11];
float* numArray5 = new float[6];
numArray1[0] = 0.3292f;
numArray1[1] = 0.254f;
numArray1[2] = 0.1021f;
numArray1[3] = -0.0073f;
numArray1[4] = -0.017f;
numArray1[5] = 0.0036f;
numArray2[0] = 3.0f / 16.0f;
numArray2[1] = 0.1704f;
numArray2[2] = 0.1266f;
numArray2[3] = 0.0735f;
numArray2[4] = 0.0298f;
numArray2[5] = 0.0059f;
numArray3[0] = 0.14272f;
numArray3[1] = 0.1217f;
numArray3[2] = 0.11401f;
numArray3[3] = 0.09303f;
numArray3[4] = 0.06583f;
numArray3[5] = 89.0f / (961.0f * 2.71828);
float* numArray6 = fUDef;
float** numArray7 = new float*[ (numrec + 1)];
for(int j=0;j<numrec+1;j++) {
numArray7[j] = new float[(numChannel + 1)];
}
float* Fil;
switch (Ftype)
{
case enumFilterType::fltWeak:
Fil = numArray1;
break;
case enumFilterType::fltNormal:
Fil = numArray2;
break;
case enumFilterType::fltStrong:
Fil = numArray3;
break;
case enumFilterType::fltUDef:
Fil = numArray6;
break;
default:
int num1 =0;// (int) Interaction.MsgBox((object) "Wrong Filter Type! Filter Set to Weak!");
Fil = numArray1;
break;
}
int num2 = numChannel;
int index1 = 1;
while (index1 <= num2)
{
int num3 = numrec;
int index2 = 1;
while (index2 <= num3)
{
if ((double) inpcurve[index2][index1] == (double) Empty)
{
numArray7[index2][index1] = Empty;
}
else
{
int index3 = -1;
int num4 = (index2 - 5);
int num5 = (index2 + 5);
int index4 = num4;
while (index4 <= num5)
{
{ ++index3; }
b[index3] = !(index4 < 1 | index4 > numrec) ? ((double) inpcurve[index4][index1] != (double) Empty ? inpcurve[index4][index1] : inpcurve[index2][index1]) : inpcurve[index2][index1];
{ ++index4; }
}
numArray7[index2][index1] = Filt(b, Fil);
}
{ ++index2; }
}
{ ++index1; }
}
for(int i=0;i<numrec+1;i++) {
for(int j=0;j<numChannel+1;j++) inpcurve[i][j] = numArray7[i][j];
}
for(int j=0;j<numrec+1;j++) delete numArray7[j];
delete numArray1;
delete numArray2;
delete numArray3;
delete numArray4;
delete b;
delete numArray5;
}
void FilterMatrixH(
float** inpcurve,
enumFilterType Ftype,
int numrec,
int numChannel,
float* fUDef,
float Empty)
{
float* numArray1 = new float[6];
float* numArray2 = new float[6];
float* numArray3 = new float[6];
float* numArray4 = new float[6];
float* b = new float[11];
float* numArray5 = new float[6];
numArray1[0] = 0.3292f;
numArray1[1] = 0.254f;
numArray1[2] = 0.1021f;
numArray1[3] = -0.0073f;
numArray1[4] = -0.017f;
numArray1[5] = 0.0036f;
numArray2[0] = 3.0f / 16.0f;
numArray2[1] = 0.1704f;
numArray2[2] = 0.1266f;
numArray2[3] = 0.0735f;
numArray2[4] = 0.0298f;
numArray2[5] = 0.0059f;
numArray3[0] = 0.14272f;
numArray3[1] = 0.1217f;
numArray3[2] = 0.11401f;
numArray3[3] = 0.09303f;
numArray3[4] = 0.06583f;
numArray3[5] = 89.0f / (961.0f * 2.71828);
float* numArray6 = fUDef;
float** numArray7 = new float*[ (numrec + 1)];
for(int j=0;j<numrec+1;j++) {
numArray7[j] = new float[(numChannel + 1)];
}
float* Fil;
switch (Ftype)
{
case enumFilterType::fltWeak:
Fil = numArray1;
break;
case enumFilterType::fltNormal:
Fil = numArray2;
break;
case enumFilterType::fltStrong:
Fil = numArray3;
break;
case enumFilterType::fltUDef:
Fil = numArray6;
break;
default:
int num1 =0;// (int) Interaction.MsgBox((object) "Wrong Filter Type! Filter Set to Weak!");
Fil = numArray1;
break;
}
int num2 = numrec;
int index1 = 1;
while (index1 <= num2)
{
int num3 = numChannel;
int index2 = 1;
while (index2 <= num3)
{
if ((double) inpcurve[index1][index2] == (double) Empty)
{
numArray7[index1][index2] = Empty;
}
else
{
int index3 = -1;
int num4 = (index2 - 5);
int num5 = (index2 + 5);
int index4 = num4;
while (index4 <= num5)
{
{ ++index3; }
b[index3] = !(index4 < 1 | index4 > numChannel) ? ((double) inpcurve[index1][index4] != (double) Empty ? inpcurve[index1][index4] : inpcurve[index1][index2]) : inpcurve[index1][index2];
{ ++index4; }
}
numArray7[index1][index2] = Filt(b, Fil);
}
{ ++index2; }
}
{ ++index1; }
}
for(int i=0;i<numrec+1;i++) {
for(int j=0;j<numChannel+1;j++) inpcurve[i][j] = numArray7[i][j];
}
for(int j=0;j<numrec+1;j++) delete numArray7[j];
delete numArray1;
delete numArray2;
delete numArray3;
delete numArray4;
delete b;
delete numArray5;
}
float Sigma(float* x,int Length, float dt, float Min)
{
float num1 = dt / 2.0f;
int num2 = 0;
float num3 = 0.0f;
float x1 = 0.0f;
float num4 = 0.0f;
float num5 = 0.0f;
int num6 = (Length - 1);
int index = 1;
while (index <= num6)
{
if ((double) x[index] > (double) Min)
{
float num7 = (float) log((double) x[index]);
float x2 = (float) index * dt - num1;
num3 += num7;
x1 += x2;
num5 += x2 * num7;
num4 += (float) pow((double) x2, 2.0);
{ ++num2; }
}
{ ++index; }
}
return !(num2 <= 0 | (double) num2 * (double) num4 - pow((double) x1, 2.0) == 0.0) ? (float) (abs(((double) num2 * (double) num5 - (double) x1 * (double) num3) / ((double) num2 * (double) num4 - pow((double) x1, 2.0))) * 4550.0) : -9999.0f;
}
float Sigma(float* x, float dt, int fromI, int toI, float Min)
{
float* x1 = new float[ (toI - fromI + 1 + 1)];
int index1 = 0;
int num1 = fromI;
int num2 = toI;
int index2 = num1;
while (index2 <= num2)
{
{ ++index1; }
x1[index1] = x[index2];
{ ++index2; }
}
float val=Sigma(x1,(toI - fromI + 1 + 1), dt, Min);
delete x1;
return val;
}
void linearSmooth7( float in[], float out[], int N, float rate[], int count ){
float *temp = new float[N];
for(int j = 0; j < N; j++)
temp[j] = in[j];
float tmpY;
for(int c = 0; c < count; c++){
int i;
if (N <= 7)
for ( i = 0; i < N; i++ ) out[i] = temp[i];
else{
out[0] = temp[0];
out[1] = temp[1];
out[2] = temp[2];
for (int i = 3; i <= N - 4; i++){
tmpY = (rate[0] * temp[i-3] + rate[1] * temp[i-2] + rate[2] *temp[i-1] + rate[3] *temp[i] + rate[4] * temp[i+1] + rate[5] * temp[i+2] + rate[6] * temp[i+3]) / 7.0;
out[i] = tmpY;
}
out[N - 3] = temp[N - 3];
out[N - 2] = temp[N - 2];
out[N - 1] = temp[N - 1];
}
for(int j = 0; j < N; j++)
temp[j] = out[j];
}
delete []temp;
}
void linearSmooth3( float in[], float out[], int N, float rate[], int count){
float *temp = new float[N];
for(int j = 0; j < N; j++)
temp[j] = in[j];
for(int c = 0; c < count; c++){
int i;
if ( N < 3 ){
for ( i = 0; i <= N - 1; i++ )
out[i] = temp[i];
}
else{
out[0] = ( 5.0 * temp[0] + 2.0 * temp[1] - temp[2] ) / 6.0;
for ( i = 1; i <= N - 2; i++ )
out[i] = ( rate[0] * temp[i - 1] + rate[1] * temp[i] + rate[2] * temp[i + 1] ) / 3.0;
out[N - 1] = ( 5.0 * temp[N - 1] + 2.0 * temp[N - 2] - temp[N - 3] ) / 6.0;
}
for(int j = 0; j < N; j++)
temp[j] = out[j];
}
delete []temp;
}
void linearSmooth5( float in[], float out[], int N, float rate[], int count ){
float *temp = new float[N];
for(int j = 0; j < N; j++)
temp[j] = in[j];
for(int c = 0; c < count; c++){
int i;
if (N <= 5)
for ( i = 0; i < N; i++ ) out[i] = temp[i];
else{
float tmpY = (3.0 * temp[0] + 2.0 * temp[1] + temp[2] - temp[4]) / 5.0;
out[0] = tmpY;
tmpY = (4.0 * temp[0] + 3.0 * temp[1] + 2 * temp[2] + temp[3]) / 10.0;
out[1] = tmpY;
for (int i = 2; i <= N - 3; i++){
tmpY = (rate[0] * temp[i - 2] + rate[1] * temp[i - 1] + rate[2] *temp[i] + rate[3] *temp[i + 1] + rate[4] * temp[i + 2]) / 5.0;
out[i] = tmpY;
}
tmpY = (4.0 * temp[N - 1] + 3.0 * temp[N - 2] + 2 * temp[N - 3] + temp[N - 4]) / 10.0;
out[N - 2] = tmpY;
tmpY = (3.0 * temp[N - 1] + 2.0 * temp[N - 2] + temp[N - 3] - temp[N - 5]) / 5.0;
out[N - 1] = tmpY;
}
for(int j = 0; j < N; j++)
temp[j] = out[j];
}
delete []temp;
}
void DataSmooth(int cal,int mode,float *in,float *out,int n,int count)
{
float rate[7];
if(cal == 0)
for(int i = 0; i < 7; i++)
rate[i] = 1;
else if(cal == 2){
if(mode == 0){rate[0] = 0.5;rate[1] = 2;rate[2] = 0.5;}
else if(mode == 1){rate[0] = 0.4;rate[1] = 0.85;rate[2] = 2.5;rate[3] = 0.85;rate[4] = 0.4;}
else if(mode == 2){rate[0] = 0.25;rate[1] = 0.5;rate[2] = 1;rate[3] = 3.5;rate[4] = 1;rate[5] = 0.5;rate[6] = 0.25;}
}
else if(cal == 3){
if(mode == 0){rate[0] = 1.25;rate[1] = 0.5;rate[2] = 1.25;}
else if(mode == 1){rate[0] = 1.3;rate[1] = 0.9;rate[2] = 0.6;rate[3] = 0.9;rate[4] = 1.3;}
else if(mode == 2){rate[0] = 1.75;rate[1] = 0.9;rate[2] = 0.6;rate[3] = 0.5;rate[4] = 0.6;rate[5] = 0.9;rate[6] = 1.75;}
}
switch(mode){
case 0:
linearSmooth3(in, out, n, rate, count);
break;
case 1:
linearSmooth5(in, out, n, rate, count);
break;
case 2:
linearSmooth7(in, out, n, rate, count);
break;
}
}
void qsort(char**datx,int left,int right,int idx,int len,int pos)
{
if(right<left+1) return;
int i,j;
float x,y;
i=left;
j=right;
x=*(float*)(&datx[(left+right)/2][0]+pos);
do
{
if(idx==0)
{
while(*(float*)(&datx[i][0]+pos)<x && i<right) i++;
while(*(float*)(&datx[j][0]+pos)>x && j>left ) j--;
}
else
{
while(*(float*)(&datx[i][0]+pos)>x && i<right) i++;
while(*(float*)(&datx[j][0]+pos)<x && j>left ) j--;
}
if(i<=j)
{
char *yy=datx[i];
datx[i]=datx[j];
datx[j]=yy;
i++;j--;
}
} while(i<=j);
if(left<j) qsort(datx,left,j,idx,len,pos);
if(i<right) qsort(datx,i,right,idx,len,pos);
}
void qsort(QStringList &datx,float *daty,int left,int right,int idx)
{
if(right<left+1) return;
int i,j;
float x,y;
i=left;
j=right;
x=datx[(left+right)/2].length();
do
{
if(idx==0)
{
while(datx[i].length()<x && i<right) i++;
while(datx[j].length()>x && j>left ) j--;
}
else
{
while(datx[i].length()>x && i<right) i++;
while(datx[j].length()<x && j>left ) j--;
}
if(i<=j)
{
QString xx=datx[i];
datx[i]=datx[j];
datx[j]=xx;
float yy=daty[i];
daty[i]=daty[j];
daty[j]=yy;
i++;j--;
}
} while(i<=j);
if(left<j) qsort(datx,daty,left,j,idx);
if(i<right) qsort(datx,daty,i,right,idx);
}
void qsort_int(char**datx,int left,int right,int idx,int len,int pos)
{
if(right<left+1) return;
int i,j;
float x,y;
i=left;
j=right;
x=*(int*)(&datx[(left+right)/2][0]+pos);
do
{
if(idx==0)
{
while(*(int*)(&datx[i][0]+pos)<x && i<right) i++;
while(*(int*)(&datx[j][0]+pos)>x && j>left ) j--;
}
else
{
while(*(int*)(&datx[i][0]+pos)>x && i<right) i++;
while(*(int*)(&datx[j][0]+pos)<x && j>left ) j--;
}
if(i<=j)
{
char *yy=datx[i];
datx[i]=datx[j];
datx[j]=yy;
i++;j--;
}
} while(i<=j);
if(left<j) qsort_int(datx,left,j,idx,len,pos);
if(i<right) qsort_int(datx,i,right,idx,len,pos);
}