/******************************************************************************* < W32IIRtoPS.c > W32の中のチャネルテーブルで選択したデータにIIRフィルタを掛けてPSファイルに出力 ********************************************************************************/ #include #include #include #include #include #define MAXCHN 60 #define MAX_IIR 16 #define DEBUG 0 #define BL 1 #define STEPKANKAKU 5.0 #define PI 3.1415926535897932384626433832795 //------ IIR 関係-------------------------- typedef unsigned char BYTE; typedef struct{ // ２次伝達関数（ｓ） double A,B,C,D,E,T; } TRNS_FNC_S; typedef struct{ // ２次伝達関数（ｚ） int sel; // 次数 0/1/2 double smp; // ｻﾝﾌﾟﾘﾝｸﾞ周波数 double a0,a1,a2,b0,b1,b2; double G; // double T; // } TRNS_FNC_Z; double ginF[MAXCHN]; // 極性反転フラグ(1=反転) int SAMP_Hz[MAXCHN]; // ｻﾝﾌﾟﾘﾝｸﾞ周波数 double ofsDt[MAXCHN]; // オフセット値 int winID0[MAXCHN],winID1[MAXCHN]; TRNS_FNC_Z IIR[MAXCHN][MAX_IIR]; // フィルタ定数バッファ //-------------------------------------- struct ChnInfo { int wch; // winﾁｬﾈﾙ番号 double sens,VpC,gin0,gin,smp; // 感度、ビット感度、ゲイン、サンプリング double abmx,max,min,avr,dif; // 絶対最大値、最大,最小、平均 (単位) double abmxV,maxV,minV,avrV,difV; // 絶対最大値、最大,最小、平均 (電圧) double avrS,avrVs; // 区間平均　2015/08/14 double dtn; // データ数 float t,t1,tw,tx; int secN,imax,imin; char stnID[8],dirc[4],unit[24]; // 観測点ID、方向コード、単位 char sTim[16],eTim[16]; // ﾁｬﾈﾙﾃｰﾌﾞﾙに書かれた開始時刻、終了時刻 time_t sTim_t,eTim_t; // int dlyT; // 開始時刻からの差(秒) char D,L; // グループ(U,N,Eの内一つ）,ライン番号0,1,2 char W[16]; char G[16]; char We[64]; int stp; // ステップ間隔 }Ch[MAXCHN]; int LUT[MAXCHN]; // Look Up table 不要 int chNmb; // 総チャンネル数 char dewin32F[MAXCHN][32]; // dewinした波形ﾌｧｲﾙ名 char filtWavF[MAXCHN][32]; // ﾌｨﾙﾀ処理した波形ﾌｧｲﾙ名 char tmpPSxyF[MAXCHN][32]; double abMax; double abmx; double GAIN; // 波形倍率 double maxScal; double dPM; // 波形極性反転(ﾁｬﾈﾙﾃｰﾌﾞﾙに記述有り) // ****** int SetChanelDt(char *inFile); // ﾁｬﾈﾙﾃｰﾌﾞﾙﾌｧｲﾙを読み保存 int MakeExIIRflt(char *filtNm); // ﾍﾙﾌﾟ機能として参考のIIRﾌｨﾙﾀﾌｧｲﾙを作成 void sjis2euc(BYTE *kj1, BYTE *kj2); void conv_sjis2euc(FILE *fi, FILE *fo); int DeWin32(char *w, int i,char *dewFile); int WaveAvrg(int ch,double *avr, int t); // 波形の平均値を得る int SetIIR2(char *FiltFile); // IIR設定ﾌｧｲﾙを読みﾌｨﾙﾀ係数をｾｯﾄする int TransFNC2_s2z(TRNS_FNC_S *s, TRNS_FNC_Z *z); // 2次伝達関数sをzに変換 int TransFNC1_s2z(TRNS_FNC_S *s, TRNS_FNC_Z *z); // 1次伝達関数sをzに変換 int Txt2IIR(char *ifn, char *ofn,int ch); // ﾃｷｽﾄ･ﾃﾞｰﾀにﾌｨﾙﾀをかける int MakePSxy(time_t stime_t,double t0,double tw); // 波形テキストファイル作成 time_t AddASCIItime(char *time0, char *time1, int sec); // ASCII時刻を加算(秒) time0 + sec => time1 time_t BCDtm2time_t(char *Tim); // 6バイトHEX時刻を time_t 時刻に変換 time_t IncSecTim(char *dTim, char *sTim); time_t SubSecTime(char *dTim, char *sTim, int sec); int LongToChr4(long ld, char *bf); int Tim12ToBCD6(char *tim, char *bcd); int SplitLin(int n,char *b,char *argv[]); // 1ﾗｲﾝの文字列をﾃﾞﾘﾐﾀ毎分割する long XtoL(char *dt); // //********************************************************************** int main(int argc, char *argv[]) { FILE *fp,*fp1,*pfp,*sfp; int ch,ch0; int pg; int S; // 1 のときはﾜｰｸﾌｧｲﾙは消さない char cmd[128],bf1[128]; char *cp; char psFname[32],psFname0[32]; char stime0[16],stime[16],etime[16]; time_t stime_t,etime_t; char gmtF[64]; double psX0,psdX,psY0,psY,psdY,psYnow,psYnow0,dPM; double t0,tw; float tD,tA; int i; time_t tst_t; // *tst* if(argc==1){ printf("--------------------------------------------------------------------------\n"); printf("win32 waveform data specified in the channel file is drawn in the ps file.\n"); printf("W32IIRtoPS tbl w32 ps flt [t1] [tw] [s]\n"); printf(" tbl : chanel tbl file(in)\n"); printf(" w32 : win32 file (in)\n"); printf(" ps : ps file (out)\n"); printf(" flt : IIR filter file(in). if 0 filter nop\n"); printf(" t1 : first time for display[sec]\n"); printf(" tw : time width for display[sec]\n"); printf(" s : if S then do not erase the work file\n"); printf(" -- 2022/03/01 ver1.2 --\n"); printf("see W32IIRtoPS ?\n"); printf("see W32IIRtoPS tbl ? flt\\nn"); return(0); } S=0; if((argc>=8)&&(*argv[7]=='S')) S=1; GAIN=1.0; maxScal=0.0; if((argc==2)&&('?'==*argv[1])){ // IIRﾌｨﾙﾀﾌｧｲﾙを表示 MakeExIIRflt(NULL); exit(0); } if(argc>=2){ if(0>=(chNmb=SetChanelDt(argv[1]))) exit(0); // ﾁｬﾈﾙﾃｰﾌﾞﾙﾌｧｲﾙを読み保存 } for(ch=0;ch=3)&&('?'==*argv[2])){ // 参考のIIRﾌｨﾙ if(argc==3){ MakeExIIRflt("tmp.flt"); exit(0); } MakeExIIRflt(argv[3]); exit(0); } if(4>argc){ printf("PS file name nothing\n"); exit(0); }else{ strncpy(psFname,argv[3],31); psFname[31]=0; strncpy(psFname0,argv[3],31); psFname0[31]=0; } for(ch=0;ch=6){ t0=(double)atof(argv[5]); }else{ t0=0.; } if(argc>=7){ tw=(double)atof(argv[6]); }else{ tw=(float)Ch[0].secN-t0; } tD=pow(10.0,(float)(int)log10(tw))/10.0; tA=tD*10.; // 開始／終了時刻が一致しているかチェック strncpy(stime0,Ch[0].sTim,14);stime0[14]=0; strncpy(etime,Ch[0].eTim,14); etime[14]=0; stime_t=AddASCIItime(stime0, stime, (int)t0); // for(ch=0;ch Ch[ch].sTim_t) Ch[ch].dlyT=(int)Ch[ch].sTim_t-(int)stime_t; else Ch[ch].dlyT=0; Ch[ch].dlyT=(int)Ch[ch].sTim_t-(int)stime_t; } // 全波形の平均値を得る for(ch=0; ch(Ch[ch].dtn=(double)WaveAvrg(ch,&(Ch[ch].avr),32767))) exit(0); } if(0!=strcmp(argv[4],"0")){ // 若しIIRﾌｧｲﾙ有りなら if(0>SetIIR2(argv[4])) exit(0); // 設定ﾌｧｲﾙを読みﾌｨﾙﾀ係数をｾｯﾄする for(ch=0; chTxt2IIR(dewin32F[ch],filtWavF[ch], ch)){ exit(0); } // ﾃｷｽﾄ･ﾃﾞｰﾀにﾌｨﾙﾀをかける strcpy(dewin32F[ch],filtWavF[ch]); // ﾌｨﾙﾀ処理した波形ﾌｧｲﾙ名をdewinした波形ﾌｧｲﾙ名に } } if(0!=MakePSxy(stime_t,t0,tw)) exit(0);// 波形テキストファイル作成 //---- GMT ---- for(ch=0;ch %s\n",psX0,psY0,psFname); psYnow=psY0; for(ch0=pg*15+1;(ch0> %s\n",0.0,-psdY,psFname); psYnow=psYnow-psdY; } ch=LUT[ch0]; fprintf(sfp,"psxy %s -P",tmpPSxyF[ch]); fprintf(sfp," -W%s",Ch[ch].W); fprintf(sfp," -R%G/%G/%G/%G",0.0,tw,Ch[ch].avrS-abmx/GAIN,Ch[ch].avrS+abmx/GAIN); fprintf(sfp," -JX"); sprintf(Ch[ch].We,"a%.8f:\"%04X[%s]\"",Ch[ch].max,Ch[ch].wch,Ch[ch].unit); // 意味不明　？ fprintf(sfp," -Ba%ef%e:\"%.2f sec\":S/we",tA,tD,tw); fprintf(sfp," -X%.1fc -Y%.1fc -O -K >> %s\n",0.0,-psdY,psFname); psYnow=psYnow-psdY; // 文字 ch0=pg*15; ch=LUT[ch0]; psY=psY0-psYnow+psdY-1.4; fprintf(sfp,"pstext -R0/70/0/10 -G%s -JX4c/2c -X-1.1c -Y%.1fc -O -K <> %s\n",Ch[ch].G,psY,psFname); fprintf(sfp,"1 1.0 7 0 8 BL %s\n",Ch[ch].unit); fprintf(sfp,"1 2.5 8 0 8 BL %s\n",Ch[ch].dirc); fprintf(sfp,"1 4.0 8 0 8 BL %s\n",Ch[ch].stnID); fprintf(sfp,"1 5.5 9 0 8 BL %04X\n",Ch[ch].wch); fprintf(sfp,"EOF\n"); psYnow=psYnow+.1; psYnow0=psYnow; for(ch0=pg*15+1;(ch0> %s\n",Ch[ch].G,-psdY,psFname); fprintf(sfp,"1 1.0 7 0 8 BL %s\n",Ch[ch].unit); fprintf(sfp,"1 2.5 8 0 8 BL %s\n",Ch[ch].dirc); fprintf(sfp,"1 4.0 8 0 8 BL %s\n",Ch[ch].stnID); fprintf(sfp,"1 5.5 9 0 8 BL %04X \n",Ch[ch].wch); fprintf(sfp,"EOF\n"); psYnow=psYnow-psdY; } ch0--; ch=LUT[ch0]; fprintf(sfp,"pstext -R0/80/0/10 -G%s -JX -X13.2c -Y%.3fc -O -K <> %s\n",Ch[ch].G,-0.10,psFname); if(Ch[ch].sens<0.0){ dPM=-1.0;}else{ dPM=1.0; } fprintf(sfp,"1 1.0 7 0 8 BL min=%+.3e:%+.3eV\n",Ch[ch].min,Ch[ch].minV*dPM); fprintf(sfp,"1 2.5 7 0 8 BL avr=%+.3e:%+.3eV\n",Ch[ch].avr,Ch[ch].avrV*dPM); fprintf(sfp,"1 4.0 7 0 8 BL max=%+.3e:%+.3eV\n",Ch[ch].max,Ch[ch].maxV*dPM); fprintf(sfp,"1 5.5 7 0 8 BL dif=%+.3e:%+.3eV\n",Ch[ch].dif,Ch[ch].difV*dPM); fprintf(sfp,"EOF\n"); psYnow=psYnow-.2; for(ch0--;ch0>=pg*15;ch0--){ if(Ch[ch].sens<0.0){ dPM=-1.0;}else{ dPM=1.0; } ch=LUT[ch0]; fprintf(sfp,"pstext -R0/80/0/10 -G%s -JX -X0.0c -Y%.3fc -O -K <> %s\n",Ch[ch].G,psdY+0.01,psFname); fprintf(sfp,"1 1.0 7 0 8 BL min=%+.3e:%+.3eV\n",Ch[ch].min,Ch[ch].minV*dPM); fprintf(sfp,"1 2.5 7 0 8 BL avr=%+.3e:%+.3eV\n",Ch[ch].avr,Ch[ch].avrV*dPM); fprintf(sfp,"1 4.0 7 0 8 BL max=%+.3e:%+.3eV\n",Ch[ch].max,Ch[ch].maxV*dPM); fprintf(sfp,"1 5.5 7 0 8 BL dif=%+.3e:%+.3eV\n",Ch[ch].dif,Ch[ch].difV*dPM); fprintf(sfp,"EOF\n"); psYnow=psYnow+psdY; } psY=psYnow0-psYnow+2.; fprintf(sfp,"pstext -R0/80/0/8 -JX25c/2c -X-13.0c -Y%.1fc -O -K <> %s\n",psY,psFname); // タイトル fprintf(sfp,"1 1 10 0 4 BL %c%c%c%c/%c%c/%c%c-%c%c:%c%c:%c%c.%03d -- %.1f sec : scale=%G [%s/div.]\n",stime[0],stime[1],stime[2],stime[3],stime[4],stime[5],stime[6],stime[7],stime[8],stime[9],stime[10],stime[11],stime[12],stime[13],(int)(round(t0*1000.)-(float)(int)t0*1000.),tw,abmx/GAIN*2.0,Ch[0].unit); if(argv[4]!=NULL){ fprintf(sfp,"1 3 10 0 8 BL filter file = %s\n",argv[4]); }else{ fprintf(sfp,"1 3 10 0 8 BL filter non \n"); } // fprintf(sfp,"1 5 12 0 8 BL \n"); fprintf(sfp,"EOF\n"); fclose(sfp); // /* *tst* euc 変換　不要？ if(0==(sfp=fopen("@@tmp.sh","r"))){ printf("can't open @@tmp.sh\n"); exit(0); } printf("kitayo psFname0=%s, gmtF%s\n",psFname0,gmtF); exit(0);// *tst* if(0==(fp1=fopen(gmtF,"w"))){ printf("can't open %s\n",gmtF); exit(0); } conv_sjis2euc(sfp, fp1); fclose(sfp); fclose(fp1); system("rm @@tmp.sh\n"); */ sprintf(cmd,"cp @@tmp.sh %s\n",gmtF); if(0==(pfp=popen(cmd,"r"))) exit(0); for(;;){ if(0==fgets(bf1,127,pfp)) break; printf("%s",bf1); } pclose(pfp); sprintf(cmd,"sh %s\n",gmtF); if(0==(pfp=popen(cmd,"r"))) exit(0); for(;;){ if(0==fgets(bf1,127,pfp)) break; printf("%s",bf1); } pclose(pfp); } // 最後のゴミ掃除 // PS作成ﾌｧｲﾙ（***.sh）を手直ししたい場合は、手動で　sh ***.sh　を実効 if(S==0){ sprintf(cmd,"rm @@*.*\n"); if(0==(pfp=popen(cmd,"r"))) exit(0); for(;;){ if(0==fgets(bf1,127,pfp)) break; } pclose(pfp); sprintf(cmd,"rm %s\n",gmtF); if(0==(pfp=popen(cmd,"r"))) exit(0); for(;;){ if(0==fgets(bf1,127,pfp)) break; printf("%s",bf1); } pclose(pfp); } exit(0); } //===================================================================== // ﾁｬﾈﾙﾃｰﾌﾞﾙﾌｧｲﾙを読み保存 2022 // 戻値=ﾁｬﾈﾙ数、引数 chTbl=ﾁｬﾈﾙﾃｰﾌﾞﾙﾌｧｲﾙ名 // [0]: 16-bit channels ID (HEX) // [1]: Recording flag // [2]: Delay time on a circuit // [3]: Station code // [4]: Motion component code // [5]: Index for monitor waveform amplitude (power of 2) // [6]: Bit size of A/D converter // [7]: Sensor sensitivity // [8]: Unit of [8] // [9]: Natural period of the seismometer (s) // [10]: Dumping constant of the sensor // [11]: Preamplificacation // [12]: LSB value // [13]: Latitude of the sensor (Positive: Northern hemisphere / Negative: Southern hemisphere) // [14]: Longitude of the sensor (Positive: Eastern hemisphere / Negative: Western hemisphere) // [15]: Altitude of the sensor // [16]: Station correction for P-wave // [17]: Station correction for S-wave // [18]: Station name (Optional) int SetChanelDt(char *chTbl) { FILE *fp; int ch,n; char bf1[256],*fctr[32]; if(0==(fp=fopen(chTbl,"r"))){ printf("filter file can't open %s\n",chTbl); return(-1); } for(ch=0,n=0;ch @tmp.tmp 2>dewin@@.tmp\n",cid,w); if(0==(fp=popen(bf,"r"))) return(-1); fgets(bf,127,fp); pclose(fp); if(0==(fp=fopen("dewin@@.tmp","r"))) return(-1); fgets(bf1,128,fp); fclose(fp); if(*bf1==0) return(-1); SplitLin(8,bf1,fctr); x=XtoL(fctr[0]); if(x!=cid) return(-1); *smp=atof(fctr[1]); // サンプリン周波数を読む if(*smp<1.0) return(-1); strncpy(stime,fctr[3],8); strncpy(stime+8,fctr[3]+9,6); *(stime+14)=0; strncpy(etime,fctr[5],8); strncpy(etime+8,fctr[5]+9,6); *(etime+14)=0; cp0=strchr(fctr[6],'['); cp1=strchr(cp0,']'); *cp1=0; *secN=atoi(cp0+1); if(0==(fp=popen("rm dewin@@.tmp\n","r"))) return(-1); for(;;){ if(0==fgets(bf,128,fp)) break; } pclose(fp); if(0==(fp=popen("rm @tmp.tmp\n","r"))) return(-1); for(;;){ if(0==fgets(bf,128,fp)) break; } pclose(fp); return(0); } //===================================================================== // w32ﾌｧｲﾙから指定したﾁｬﾈﾙの波形ﾃﾞｰﾀﾌｧｲﾙを作る // 引数： winﾌｧｲﾙ名、n ,dewFile=波形ﾃﾞｰﾀﾌｧｲﾙ名　2022/02/18 int DeWin32(char *w32, int ch,char *dewFile) { FILE *fp; int i,n,x,st,et; char bf[128],*fctr[10],*cp0,*cp1; char tbf[16]; sprintf(bf,"dewin_32 -c %04X %s > %s 2>@@dewinMsg.txt\n",Ch[ch].wch,w32,dewFile); if(0==(fp=popen(bf,"r"))) return(-1); for(i=0; ; i++){ if(0==fgets(bf,127,fp)) break; } pclose(fp); // dewin_32 ﾒｯｾｰｼﾞを読む if(0==(fp=fopen("@@dewinMsg.txt","r"))){ printf("cnan't open dewin_32 message as %04X\n",Ch[ch].wch); fclose(fp); return(-1);} if(0==fgets(bf,127,fp)){ printf("dewin_32 message nothing as %04X\n",Ch[ch].wch); fclose(fp); return(-1);} fclose(fp); SplitLin(8,bf,fctr); // サンプリン周波数を読む Ch[ch].smp=(double)atoi(fctr[1]); if(Ch[ch].smp<1.0){ printf("smpling Hz 0 as %04X\n",Ch[ch].wch); return(-1); } // 開始／終了時刻を保存 strncpy(Ch[ch].sTim,fctr[3],8); strncpy(Ch[ch].sTim+8,fctr[3]+9,6); *(Ch[ch].sTim+14)=0; strncpy(Ch[ch].eTim,fctr[5],8); strncpy(Ch[ch].eTim+8,fctr[5]+9,6); *(Ch[ch].eTim+14)=0; Ch[ch].sTim_t=AddASCIItime(Ch[ch].sTim,tbf, 0); Ch[ch].eTim_t=AddASCIItime(Ch[ch].eTim,tbf, 0); // 記録時間を保存 cp0=strchr(fctr[6],'['); cp1=strchr(cp0,']'); *cp1=0; Ch[ch].secN=atoi(cp0+1); // ﾒｯｾｰｼﾌｧｲﾙを削除 return(0); } //================================================= // 波形の平均値を得る int WaveAvrg(int ch,double *avr, int t) { FILE *fp; char bf[64]; double tdt,dt; int i; if(0==(fp=fopen(dewin32F[ch],"r"))){ printf("can't open %s,dewin(%04X)\n",dewin32F[ch],Ch[ch].wch); fclose(fp);return(-1); } for(i=0,tdt=0.0;;i++){ if(0==fgets(bf,63,fp)){ break; } if((t!=0)&&(i>=t)) break; tdt=tdt+(double)atoi(bf); } fclose(fp); if(t!=0){ *avr=tdt/(double)t; } else{ *avr=tdt/(double)i; } return(i); } //================================================================= // IIR設定ﾌｧｲﾙを読みﾌｨﾙﾀ係数をｾｯﾄする int SetIIR2(char *FiltFile) { FILE *fp; char bf[256],*fctr[16]; double f,h,fh,fl; double smpHz,avl; // int chnID; int ch,i,e,n,t,sel; TRNS_FNC_S s; char *c; fl=0.0001; if(0==(fp=fopen(FiltFile,"r"))){ printf("can't open %s\n",FiltFile); exit(0); } do{ if(0==fgets(bf,255,fp)){ printf("error : %s\n",FiltFile); exit(0); } n=SplitLin(6,bf,fctr); if(0==strcmp(fctr[0],"WIN")) break; }while(*bf=='#'); for(e=0,ch=0;ch%s,%d\n",fctr[0],ch+1); e=-1; break; } SAMP_Hz[ch]=(int)Ch[ch].smp;// ｻﾝﾌﾟﾘﾝｸﾞ周波数 smpHz=Ch[ch].smp; s.T=1.0/smpHz; if((n>=4)&&(*fctr[3]!=0)) ginF[ch]=atof(fctr[3]); // GAIN if(ginF[ch]==0.0) ginF[ch]=1.0; Ch[ch].gin=Ch[ch].gin*ginF[ch]; t=0; if((n>=6)&&(*fctr[4]=='T')) t=atoi(fctr[5]); if(0>WaveAvrg(ch,&avl,t)) return(-1); // 波形の平均値を得る 2022 do{ if((n>=6)&&(*fctr[4]=='T')){ ofsDt[ch]=avl; break; } if((n>=5)&&(*fctr[4]=='A')){ ofsDt[ch]=avl; break; } if((n>=6)&&(*fctr[4]=='S')){ ofsDt[ch]=atof(fctr[5]); break; } }while(0); fh=smpHz; for(i=0; i=3) IIR[ch][i].G=IIR[ch][i].G*atof(fctr[2]); strcat(Ch[ch].unit,"/s"); break; case 2: printf("XXerror : %s,%s,%d,%d\n",FiltFile,bf,ch,i+1); e=-1; strcat(Ch[ch].unit,"/s/s"); break; default: printf("error : %s,%s,sel=%d,ch=%d,%d\n",FiltFile,bf,sel,ch,i+1); e=-1; break; } break; } if(0==strcmp(fctr[0],"DIF")){ // ** 微分 ** switch(sel){ case 0: break; case 1: s.A=2.*PI*fh; s.B=2.*PI*fh; s.C=0.0; TransFNC1_s2z(&s,IIR[ch]+i); if(n>=3) IIR[ch][i].G=IIR[ch][i].G*atof(fctr[2]); strcat(Ch[ch].unit,"*s"); break; case 2: printf("error : %s,%s,%d,%d\n",FiltFile,bf,ch,i+1); e=-1; strcat(Ch[ch].unit,"*s*s"); break; default: printf("error : %s,%s,%d,%d\n",FiltFile,bf,ch,i+1); e=-1; break; } break; } printf("error. bad filter: %s,%s\n",FiltFile,bf); e=-1; break; }while(0); if(e!=0){ break; } }// end of for i if(e!=0){ break; } if(0==strcmp(fctr[0],"WIN")){ continue; } if(0==fgets(bf,255,fp)){ break; } n=SplitLin(6,bf,fctr); } // end of for ch fclose(fp); return(e); } //================================================================= // // C･s^2 + D･s + E // Hs = --------------------- // s^2 + A･s + B // //a =2/T , G = 1/(a^2 + A*a + B), b=C*a^2 // // (b + D･a + E) + 2･(E - b)･z^-1 + (b - D･a + E)･z^-2 // Hz =G･-------------------------------------------------------- // 1 + 2･G･(B - a^2)･z^-1 + G･(a^2 - A･a + B)･z^-2 // // a2･z^-2 + a1･z^-1 + a0 // Hz =G･-------------------------- // b2･z^-2 + b1･z^-1 + b0 // int TransFNC2_s2z(TRNS_FNC_S *s, TRNS_FNC_Z *z) // 2次伝達関数sをzに変換 { double a,b,G; a = 2.0/s->T; b = s->C*(a*a); G = 1.0/((a*a) + s->A*a + s->B); z->a0 = b + s->D*a + s->E; z->a1 = 2.0*(s->E - b); z->a2 = b - s->D*a + s->E; z->b0 = 1.0; z->b1 = 2.0*G*(s->B - a*a); z->b2 = G*((a*a) - s->A*a + s->B); z->G = G; z->T = s->T; return(0); } //================================================================= // Bs + C // Hs = ----------- // s + A // // 1 (B*a + C)+(C - B*a)*z^-1 // Hz = ----------*---------------------------- // A + a A - a // 1 + --------*z^-1 // A + a//5 // // a0 + a1*z^-1 // Hz = G*------------------- // b0 + b1*z^-1 // a = 2/T // int TransFNC1_s2z(TRNS_FNC_S *s, TRNS_FNC_Z *z) // 1次伝達関数sをzに変換 { double a,G; a = 2.0/s->T; // 2/T G = 1.0/(s->A + a); // 1/(A + a) z->a0= s->B * a + s->C; // (B*a + C) z->a1= s->C - s->B*a; // (C - B*a) z->b0= 1.0; // 1 z->b1= (s->A - a)/(s->A + a); // (A - a)/(A + a) z->G = G; z->T = s->T; return(0); } //================================================================= int Txt2IIR(char *ifn, char *ofn,int ch) // ﾃｷｽﾄ･ﾃﾞｰﾀにﾌｨﾙﾀをかける { FILE *ifp,*ofp; TRNS_FNC_Z *z; int i,j,dt0; char bf[128]; double dt,Da[MAX_IIR][2],Db[MAX_IIR][2],X,Y,Dt; X=Y=dt=0.0; for(i=0;isel){ case 0 : break; case 1: // 一次 X=(z->a1*Da[j][0] + z->a0*dt); Y=(X - z->b1*Db[j][0]); Dt=z->G*Y; Da[j][0]=dt; Db[j][0]=Y; dt=Dt; break; case 2: // 二次 case -2: X=z->a2*Da[j][1] + z->a1*Da[j][0] + z->a0*dt; Y=(X - z->b2*Db[j][1] - z->b1*Db[j][0]); Dt=z->G*Y; Da[j][1]=Da[j][0]; Da[j][0]=dt; Db[j][1]=Db[j][0]; Db[j][0]=Y; dt=Dt; break; default: break; } } fprintf(ofp,"%d\n",(int)dt); } fclose(ifp); fclose(ofp); return(0); } //===================================================================== // 時刻付き波形テキストファイル作成 int MakePSxy(time_t stime_t, double t1, double tw) { FILE *fp,*fp1; char bf1[128]; int i,n,ch,imax,imin,idt; // time_t ts_t; // double t,t1,t2,td,dt,tdt,tdd; double t,t2,td,dt,tdt,tdd; double t1s,t2s,tws,tdtS,tnS; // 区間トータル2015/08/15 ○ double psX0,psdX,psY0,psY,psdY,psYnow,psYnow0,dPM; float stpT0,stpT1; int stpMx,stpMi,dStp,stp; int maxCh; int nop; t1s=t1; tws=tw; // 20145/08/15 ○ abmx=0.0; // 2022/02/23 for(ch=0;ch=0.0){ dPM=1.0; }else{ dPM=-1.0; } if(0==(fp=fopen(dewin32F[ch],"r"))){ printf("cant't open %s,\n",dewin32F[ch]); exit(0); } if(0==(fp1=fopen(tmpPSxyF[ch],"w"))){ printf("cant't open %s,\n",tmpPSxyF[ch]); exit(0); } if(Ch[ch].smp*(double)tw<(115.*STEPKANKAKU)){ // サンプル数小 imin=0x7FFFFFFF; imax=0-imin; td=1.0/Ch[ch].smp; for(t=0.0,n=0,tdt=0,i=0;;i++){ if(0==fgets(bf1,127,fp)) break; if(((int)(t1*Ch[ch].smp)<=i)&&(i<(int)((t1+tw)*Ch[ch].smp))){ idt=atoi(bf1)*(int)dPM; tdt=tdt+(float)idt; if(idt>imax){ imax=idt; } if(idt=Ch[ch].avrV-Ch[ch].minV){ Ch[ch].abmxV=Ch[ch].maxV-Ch[ch].avrV; }else{ Ch[ch].abmxV=Ch[ch].avrV-Ch[ch].minV; } Ch[ch].abmx=Ch[ch].abmxV/Ch[ch].sens*dPM; if(Ch[ch].abmx>abmx) abmx=Ch[ch].abmx; fclose(fp1); fclose(fp); }else{ //ｻﾝﾌﾟﾙ数大の為GMT描画時間大、波形を間引き(重ねた線)する imin=0x7FFFFFFF; imax=0-imin; td=1.0/Ch[ch].smp; Ch[ch].stp=(int)(Ch[ch].smp*tw/(115. * STEPKANKAKU * 2.)); // ステップ間隔 stp=0; stpMi=0x7FFFFFFF; stpMx=-stpMi-1; // 区間最小最大＝初期値 for(t=0.0,n=0,tdt=0,i=0;;i++){ if(0==fgets(bf1,127,fp)) break; if(((int)(t1*Ch[ch].smp)<=i)&&(i<(int)((t1+tw)*Ch[ch].smp))){ idt=atoi(bf1)*(int)dPM; tdt=tdt+(float)idt; if(idt>imax){ imax=idt; } if(idtstpMx){ stpMx=idt;} if(idtCh[ch].stp){ dt=(double)stpMi*Ch[ch].VpC/Ch[ch].sens/Ch[ch].gin*dPM; // fprintf(fp1,"%f %G\n",stpT0,dt); fprintf(fp1,"%f,%G\n",stpT0,dt); dt=(double)stpMx*Ch[ch].VpC/Ch[ch].sens/Ch[ch].gin*dPM; // fprintf(fp1,"%f %G\n",stpT1,dt); fprintf(fp1,"%f,%G\n",stpT1,dt); stp=0; stpMi=0x7FFFFFFF; stpMx=-stpMi-1; // 区間最小最大＝初期値 } t=t+td; n++; } if(((int)(t1s*Ch[ch].smp)<=i)&&(i<(int)((t1s+tws)*Ch[ch].smp))){ idt=atoi(bf1)*(int)dPM; tdtS=tdtS+(float)idt; tnS=tnS+1.0; } } Ch[ch].dtn=(double)n; Ch[ch].imax=imax; Ch[ch].imin=imin; Ch[ch].avrV =tdt/Ch[ch].dtn * Ch[ch].VpC/Ch[ch].gin; Ch[ch].avrVs=tdtS/tnS * Ch[ch].VpC/Ch[ch].gin; Ch[ch].maxV=(double)imax*Ch[ch].VpC/Ch[ch].gin; Ch[ch].minV=(double)imin*Ch[ch].VpC/Ch[ch].gin; Ch[ch].difV=Ch[ch].maxV-Ch[ch].minV; Ch[ch].avr=Ch[ch].avrV/Ch[ch].sens*dPM; Ch[ch].avrS=Ch[ch].avrVs/Ch[ch].sens*dPM; Ch[ch].max=Ch[ch].maxV/Ch[ch].sens*dPM; Ch[ch].min=Ch[ch].minV/Ch[ch].sens*dPM; Ch[ch].dif=Ch[ch].difV/Ch[ch].sens*dPM; if(Ch[ch].maxV-Ch[ch].avrV>=Ch[ch].avrV-Ch[ch].minV){ Ch[ch].abmxV=Ch[ch].maxV-Ch[ch].avrV; }else{ Ch[ch].abmxV=Ch[ch].avrV-Ch[ch].minV; } Ch[ch].abmx=Ch[ch].abmxV/Ch[ch].sens*dPM; if(Ch[ch].abmx>abmx){ abmx=Ch[ch].abmx; maxCh=ch; } fclose(fp1); fclose(fp); } } if(abmx==0.0){ abmx=1.0; } return(0); } //===================================================================== // ASCII時刻を加算(秒) time0 + sec => time1 time_t AddASCIItime(char *time0, char *time1, int sec) { struct tm Tm,*tM; time_t timer; char tbf[6]; tbf[2]=0; tbf[4]=0; Tm.tm_isdst= -1; // 夏時間無効 ※重要 strncpy(tbf,time0+12,2); Tm.tm_sec=atoi(tbf); strncpy(tbf,time0+10,2); Tm.tm_min=atoi(tbf); strncpy(tbf,time0+8,2); Tm.tm_hour=atoi(tbf); strncpy(tbf,time0+6,2); Tm.tm_mday=atoi(tbf); strncpy(tbf,time0+4,2); Tm.tm_mon=atoi(tbf)-1; strncpy(tbf,time0+0,4); Tm.tm_year=atoi(tbf)-1900; timer=mktime(&Tm); timer=timer+(time_t)sec; tM = localtime(&timer); sprintf(time1,"%04d%02d%02d%02d%02d%02d",tM->tm_year+1900,tM->tm_mon+1,tM->tm_mday,tM->tm_hour,tM->tm_min,tM->tm_sec); return(timer); } //===================================================================== // 6バイトHEX時刻を time_t 時刻に変換 time_t BCDtm2time_t(char *Tim) { time_t st; struct tm sT; sT.tm_isdst = -1; sT.tm_sec = (*(Tim+10) & 0x0F)*10+(*(Tim+11) & 0x0F); sT.tm_min = (*(Tim+8) & 0x0F)*10+(*(Tim+9) & 0x0F); sT.tm_hour= (*(Tim+6) & 0x0F)*10+(*(Tim+7) & 0x0F); sT.tm_mday= (*(Tim+4) & 0x0F)*10+(*(Tim+5) & 0x0F); sT.tm_mon = (*(Tim+2) & 0x0F)*10+(*(Tim+3) & 0x0F) - 1; sT.tm_year= (*(Tim+0) & 0x0F)*10+(*(Tim+1) & 0x0F); if(sT.tm_year<50){ sT.tm_year =sT.tm_year+100; } st = mktime(&sT); return(st); } //===================================================================== // yyyymmddhhmm+1 time_t IncSecTim(char *dTim, char *sTim) { time_t st,dt; struct tm sT,*dT; // int ye,mo,dy,ho,mi,cy; sT.tm_isdst = -1; sT.tm_sec = (*(sTim+10) & 0x0F)*10+(*(sTim+11) & 0x0F); sT.tm_min = (*(sTim+8) & 0x0F)*10+(*(sTim+9) & 0x0F); sT.tm_hour= (*(sTim+6) & 0x0F)*10+(*(sTim+7) & 0x0F); sT.tm_mday= (*(sTim+4) & 0x0F)*10+(*(sTim+5) & 0x0F); sT.tm_mon = (*(sTim+2) & 0x0F)*10+(*(sTim+3) & 0x0F) - 1; sT.tm_year=(*(sTim+0) & 0x0F)*10+(*(sTim+1) & 0x0F); if(sT.tm_year<50){ sT.tm_year =sT.tm_year+100; } st = mktime(&sT); dt = st +(time_t)1; dT = localtime(&dt); sprintf(dTim,"%02d%02d%02d%02d%02d%02d", dT->tm_year-100,dT->tm_mon+1,dT->tm_mday,dT->tm_hour,dT->tm_min,dT->tm_sec); return(dt); } //------------------------------------------------------------------------------- time_t SubSecTime(char *dTim, char *sTim, int sec) { time_t st,dt; struct tm sT,*dT; sT.tm_isdst = -1; sT.tm_sec = (*(sTim+10) & 0x0F)*10+(*(sTim+11) & 0x0F); sT.tm_min = (*(sTim+8) & 0x0F)*10+(*(sTim+9) & 0x0F); sT.tm_hour= (*(sTim+6) & 0x0F)*10+(*(sTim+7) & 0x0F); sT.tm_mday= (*(sTim+4) & 0x0F)*10+(*(sTim+5) & 0x0F); sT.tm_mon = (*(sTim+2) & 0x0F)*10+(*(sTim+3) & 0x0F) - 1; sT.tm_year=(*(sTim+0) & 0x0F)*10+(*(sTim+1) & 0x0F); if(sT.tm_year<50){ sT.tm_year =sT.tm_year+100; } st = mktime(&sT); dt = st -(time_t)sec; dT = localtime(&dt); sprintf(dTim,"%02d%02d%02d%02d%02d%02d", dT->tm_year-100,dT->tm_mon+1,dT->tm_mday,dT->tm_hour,dT->tm_min,dT->tm_sec); return(dt); } //------------------------------------------------------------------------------- int LongToChr4(long ld, char *bf) { bf[0]=ld >> 24; bf[1]=(ld & 0x00FFFFFF) >> 16; bf[2]=(ld & 0x0000FFFF) >> 8; bf[3]=(ld & 0x000000FF); return(0); } //------------------------------------------------------------------------------- int Tim12ToBCD6(char *tim, char *bcd) { char *spo,*dpo; int i; spo=tim; dpo=bcd; for(i=0;i<6;i++){ *dpo=(*spo & 0x0F)*16 + (*(spo+1) & 0x0F); spo=spo+2; dpo++; } return(0); } //===================================================================== // 1ﾗｲﾝの文字列をﾃﾞﾘﾐﾀ毎分割する 03/06/19 int SplitLin(int max_n,char *bf,char *argv[]) { int n; char *cp; n=0; cp=bf; argv[0]=cp; while(*cp!=0){ switch(*cp){ case 0: return(n); case 0x0d: *cp=0; return(n); case 0x0a: *cp=0; return(n); case ' ': // 先頭がｽﾍﾟｰｽなら while(*cp==' '){ // ｽﾍﾟｰｽが続くかぎり *cp=0; // 0を書込 cp++; // ﾎﾟｲﾝﾀを++ if(*cp==0) return(n); } break; case ',': // ｶﾝﾏなら while(*cp==','){ *cp=0; cp++; // 0を書く、ﾎﾟｲﾝﾀを++ if(*cp!=',') break; argv[n]=cp; n++; // ｱﾄﾞﾚｽを保存、ﾎﾟｲﾝﾀを++ if(n>=max_n) return(n); if(*cp==0) return(n); } break; default: // それ以外なら argv[n]=cp; n++; // 先頭ｱﾄﾞﾚｽを保存、ﾎﾟｲﾝﾀを++ if(n>=max_n) return(n); while((*cp!=' ')&&(*cp!=',')){ if((*cp==0x0a)||(*cp==0x0d)){ *cp=0; return(n); } if(*cp==0) return(n); cp++; }; break; } } return(n); } //------------------------------------------------------------------------------- long XtoL(char *dt) { int i,j,l,n; long a,c; l=strlen(dt); a=0; n=0; for(i=l-1;i>=0;i--){ c=*(dt+i); if(c<'0') return(-1); if(c<='9'){ c =c & 0x0f; }else{ c =c & 0x0f; c =c + 9; } for(j=0;j 0x7E) *kj2 += 0x60; else *kj2 += 0x61; } else { if (*kj1 < 0x3F) *kj1 -= 0x60; else *kj1 += 0x20; *kj2 += 0x02; } } //------------------------------------------------------------------------------- void conv_sjis2euc(FILE *fi, FILE *fo) { int n; BYTE c1, c2 = 0; while ((n = getc(fi)) != EOF) { c1 = (BYTE)n; if (c1 < 0x80) { fputc(c1, fo); continue; } if (0xA1 <= c1 && c1 <= 0xDF) { fputc(0x8E, fo); fputc(c1, fo); continue; } if ((n = getc(fi)) == EOF) break; c2 = (BYTE)n; sjis2euc(&c1, &c2); fputc(c1, fo); fputc(c2, fo); } }