00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023 #include <stdlib.h>
00024 #include <stdio.h>
00025 #include <time.h>
00026 #include <ctype.h>
00027 #include <math.h>
00028 #include <errno.h>
00029 #include <unistd.h>
00030 #include <string.h>
00031 #include <sys/types.h>
00032 #include <sys/stat.h>
00033 #include <fcntl.h>
00034
00035 #include "common.h"
00036
00037 int32_t msf[9][3] = {
00038 { 0, 0, 0 },
00039 { 0xffffff, 1, 8 },
00040 { 0x000fff, 12, 20 },
00041 { 0x0000ff, 8, 24 },
00042 { 0x00003f, 6, 26 },
00043 { 0, 0, 0 },
00044 { 0x00000f, 4, 28 },
00045 { 0, 0, 0 },
00046 { 0x000007, 3, 29 }
00047 };
00048
00049 int32_t lastData[MAX_CHANNEL]={
00050 0, 0, 0, 0, 0, 0, 0, 0,
00051 0, 0, 0, 0, 0, 0, 0, 0,
00052 0, 0, 0, 0, 0, 0, 0, 0,
00053 0, 0, 0, 0, 0, 0, 0, 0,
00054 0, 0, 0, 0, 0, 0, 0, 0
00055 };
00056
00057 int nTitan1=0,nTitan2=0;
00058 int lastFormat=-1;
00059
00060 struct Titan2TimeType lastFrameTime = {0, 0};
00061
00062 int64_t lastSegment=0;
00063
00064
00065
00066
00067 void SetAuxChannel1(int ch, int32_t val)
00068 {
00069 HWConfig.data[ch][HWConfig.srObs[ch].npts]=val;
00070 HWConfig.srObs[ch].npts++;
00071 HWConfig.srExp[ch].npts++,
00072 HWConfig.srExp[ch].base=HWConfig.srExp[0].base;
00073 HWConfig.srExp[ch].div= HWConfig.srExp[0].div*128;
00074 }
00075
00076 void ProcessAux1Frame(const unsigned char *frame)
00077 {
00078 int pair,ch;
00079 int32_t V0;
00080 int32_t V1;
00081
00082 if (auxiliaryMask==0LL)
00083 {
00084 return;
00085 }
00086
00087 pair=frame[9]&0x7;
00088
00089 V0=frame2i24(frame,0);
00090 V1=frame2i24(frame,3);
00091
00092 ch=MAX_OSIRIS_CHANNEL+pair*2+0;
00093 PrintDebug("AUX: pair=%d ch=%d V0=%d V1=%d\n",pair,ch,V0,V1);
00094 SetAuxChannel1(ch+0, V0);
00095 SetAuxChannel1(ch+1, V1);
00096 }
00097
00098
00099
00100
00101 void SetAuxChannel2(int ch, int32_t val)
00102 {
00103 HWConfig.data[ch][0]=val;
00104 HWConfig.srObs[ch].npts=1;
00105 HWConfig.srExp[ch].npts=1;
00106 HWConfig.srExp[ch].base=1;
00107 HWConfig.srExp[ch].div=1;
00108 }
00109
00111 void FlushDataFiles(void)
00112 {
00113 int i;
00114 for (i=0; i<MAX_CHANNEL; i++)
00115 {
00116 CloseMSeedFile(i);
00117 CloseSacFile(i);
00118 CloseSegyFile(i);
00119 CloseSeg2File(i);
00120 CloseAscFile(i);
00121 CloseBinFile(i);
00122 CloseWavFile(i);
00123 }
00124 CloseFullSeg2File();
00125 CloseGifFile();
00126 CloseTitFile();
00127 CloseSisFile();
00128
00129 lastFrameTime.usecond=0;
00130 }
00131
00133 void FlushNonDataFiles(void)
00134 {
00135 CloseTimeFile(oldHWConfig.fieldIdent);
00136 CloseCacheFile(0);
00137 }
00138
00147 int CheckContinuity(void)
00148 {
00149 int i,it=-1,ret;
00150 double eDeltaT;
00151 double oDeltaT;
00152 double threshold;
00153 struct Titan2TimeType diffTime = { 0, 0 };
00154 ret=0;
00155 if (ignoreDiscontinuity)
00156 {
00157 threshold=120.;
00158 }
00159 else
00160 {
00161 threshold=1./900.;
00162 }
00163
00164
00165 for (i=0; i<MAX_CHANNEL; i++)
00166 {
00167
00168 if (
00169 (HWConfig.srExp[i].base<1)||
00170 (HWConfig.srExp[i].div<1)||
00171 (HWConfig.srExp[i].npts<1)
00172 )
00173 {
00174 continue;
00175 }
00176
00177 if (HWConfig.srExp[i].base/HWConfig.srExp[i].div<1)
00178 {
00179 continue;
00180 }
00181
00182 if (it==-1)
00183 {
00184 it=i;
00185 }
00186
00187 if (HWConfig.srObs[i].npts!=HWConfig.srExp[i].npts)
00188 {
00189 PrintError("ch%d, obs/exp=%d/%d (%s @ %d)%s\n",
00190 i,
00191 HWConfig.srObs[i].npts,
00192 HWConfig.srExp[i].npts,
00193 currentFile,SafeFtell(titan2File),
00194 (HWConfig.srObs[i].npts<HWConfig.srExp[i].npts)?"":
00195 " killed");
00196 if (HWConfig.srObs[i].npts>HWConfig.srExp[i].npts)
00197 {
00198 PrintDebug("killing the frame: it is too long\n");
00199
00200 HWConfig.srObs[i].npts=0;
00201 HWConfig.srExp[i].npts=0;
00202 lastData[i]=0;
00203 }
00204
00205 if (extractingData)
00206 {
00207 ret=-2;
00208 }
00209 else
00210 {
00211 HWConfig.srObs[i].npts=0;
00212 HWConfig.srExp[i].npts=0;
00213 lastData[i]=0;
00214 }
00215 }
00216 }
00217 if (it==-1)
00218 {
00219 return(0);
00220 }
00221
00222
00223 if ((nTitan1+nTitan2)==0)
00224 {
00225 return(ret);
00226 }
00227
00228
00229 if ((lastFrameTime.usecond>=(1LL<<32))&&(currentTime.usecond>(1LL<<32)))
00230 {
00231 eDeltaT=DiffTimeT2T1(&lastFrameTime,¤tTime,&diffTime);
00232 lastFrameTime.usecond=currentTime.usecond;
00233 if ((it!=-1)&&(HWConfig.srExp[it].base>0))
00234 {
00235 oDeltaT=(double)(HWConfig.srObs[it].npts*HWConfig.srExp[it].div)/
00236 (double)HWConfig.srExp[it].base;
00237 }
00238 else
00239 {
00240 oDeltaT=0.;
00241 }
00242
00243 if ((HWConfig.formatVersion==1)&&(extractingData)&&(it<6))
00244 {
00245 if (fabs((oDeltaT-eDeltaT)-
00246 (float)HWConfig.srExp[it].div/(float)HWConfig.srExp[it].base)<0.001)
00247 {
00248 PrintWarning("one sample corruption detected (fixed): %s @ %d\n\t%s\n",
00249 currentFile, SafeFtell(titan2File), StrTime(¤tTime));
00250 for (i=0; i<MAX_OSIRIS_CHANNEL; i++)
00251 {
00252 int _i;
00253 int _imax=1,vtest=0,vmax=0;
00254 if (HWConfig.srExp[i].npts<1) continue;
00255 for (_i=1; _i<127; _i++)
00256 {
00257 vtest=abs(HWConfig.data[i][_i]-HWConfig.data[i][_i-1])/(1+abs(abs(HWConfig.data[i][_i-1])-abs(HWConfig.data[i][_i+1])));
00258 if (vtest>vmax)
00259 {
00260 vmax=vtest;
00261 _imax=_i;
00262 }
00263
00264
00265
00266
00267
00268
00269 }
00270
00271
00272
00273
00274
00275
00276
00277
00278
00279 HWConfig.data[i][_imax]=(HWConfig.data[i][_imax-1]+HWConfig.data[i][_imax+1])/2;
00280 }
00281
00282 eDeltaT+=(float)HWConfig.srExp[it].div/(float)HWConfig.srExp[it].base;
00283 }
00284
00285 }
00286
00287 if (fabs(eDeltaT-oDeltaT)>threshold)
00288 {
00289
00290 ChangeStreamConfig();
00291 if (seeking==0)
00292 {
00293 PrintWarning("data %s: %s @ %d\n\t%s "
00294 "%+.4fs (TimeFrame:%.3fs, Samples:%.3fs)\n",
00295 (eDeltaT>oDeltaT)?"lag":"overlap",
00296 currentFile,
00297 SafeFtell(titan2File),
00298 StrTime(¤tTime),
00299 eDeltaT-oDeltaT,
00300 eDeltaT,oDeltaT);
00301 LogPrint("DISCONTINUITY(s) %f\n",eDeltaT-oDeltaT);
00302 }
00303 if (ret>=0)
00304 {
00305 ret=-1;
00306 }
00307
00308 if (HWConfig.srExp[it].npts==128)
00309 {
00310 int ed;
00311 double f;
00312 f=(double)HWConfig.srExp[it].base/(double)HWConfig.srExp[it].div;
00313 ed=(int)((100000000.1*(eDeltaT-oDeltaT)*f)/100000.0);
00314 PrintWarning("%f %d\n",
00315 100000000.1*(eDeltaT-oDeltaT)*f,
00316 ed);
00317 if (ed%(HWConfig.srExp[it].npts*1000)==0)
00318 {
00319 PrintWarning("POSSIBLE SUPERFRAME JUMP: %d frames\n",
00320 ed/HWConfig.srExp[it].npts/1000);
00321 LogPrint("SuperFrameSkip %d\n",
00322 ed/HWConfig.srExp[it].npts/1000);
00323 }
00324 }
00325 else
00326 {
00327 PrintWarning("npts=%d\n",HWConfig.srExp[it].npts);
00328 }
00329 }
00330 else
00331 {
00332
00333 if (fabs(eDeltaT-oDeltaT)>1./500.)
00334 {
00335 if (PrintMessages())
00336 {
00337 PrintWarning("time uncertainity: %s @ %d\n\t%s %+.4fs\n",
00338 currentFile,
00339 SafeFtell(titan2File),
00340 StrTime(¤tTime),
00341 eDeltaT-oDeltaT);
00342 }
00343 }
00344 }
00345 }
00346 lastFrameTime.usecond=currentTime.usecond;
00347 return(ret);
00348 }
00349
00354 void EndSuperFrame(void)
00355 {
00356 int i,continuity=0;
00357 int exitNow=0;
00358 int score,format;
00359
00360
00361
00362 continuity=CheckContinuity();
00363
00364
00365 if (timeSegment>0)
00366
00367 {
00368 int64_t seg;
00369 beginWindow=0;
00370
00371 seg=correctedTime.usecond/timeSegment;
00372 seg>>=32;
00373 if (seg!=lastSegment)
00374 {
00375 endWindow=1;
00376 lastSegment=seg;
00377 }
00378 else
00379 {
00380 endWindow=0;
00381 }
00382
00383 }
00384
00385 exitNow=0;
00386
00387
00388 if (endWindow>0)
00389 {
00390 if (scanAll==0)
00391 {
00392 if (dataExtracted)
00393 {
00394 if ((correctedTime.usecond>>32)>endWindow)
00395 {
00396 exitNow=1;
00397 }
00398 }
00399 }
00400 }
00401
00402 if (onlyTag[0])
00403 {
00404 if ((scanAll==0)&&(strcmp(currentTag,onlyTag)!=0))
00405 {
00406 exitNow=1;
00407 }
00408 }
00409
00410
00411 if (
00412 (endWindow>0)&&
00413 (beginWindow>0)&&
00414 ((correctedTime.usecond>>32)>endWindow)&&
00415 ((correctedTime.usecond>>32)>beginWindow)
00416 )
00417 {
00418 PrintDebug("endWindow=%d beginWindow=%d t=%d (%s)\n",
00419 endWindow,beginWindow,
00420 (int)((correctedTime.usecond>>32)),
00421 StrTime(&correctedTime));
00422 if (scanAll==0)
00423 {
00424 exitNow=1;
00425 }
00426 }
00427
00428 if (exitNow)
00429 {
00430
00431 if (extractingData)
00432 {
00433 PrintLog("end of window detected\n");
00434 }
00435
00436 extractingData=0;
00437 FlushDataFiles();
00438
00439 if (timeSegment>0)
00440
00441 {
00442
00443 NextWindow();
00444
00445 }
00446 else
00447 {
00448
00449 if ((timeList>0)||((beginWindow>0)&&(endWindow>0)))
00450
00451 {
00452
00453 CloseInput();
00454
00455 }
00456
00457 }
00458 }
00459
00460
00461
00462
00463
00464
00465
00466
00467
00468
00469
00470
00471
00472
00473
00474
00475
00476
00477
00478
00479
00480
00481
00482
00483
00484
00485
00486
00487
00488
00489
00490
00491
00492
00493
00494
00495
00496
00497
00498
00499
00500
00501
00502
00503
00504
00505
00506
00507
00508
00509
00510
00511
00512
00513
00514
00515
00516
00517
00518
00519
00520
00521
00522
00523
00524
00525
00526
00527
00528
00529
00530
00531
00532 if ((continuity<0)&&(extractTit==0))
00533 {
00534 PrintDebug("closing data files\n");
00535 CachePrintCST();
00536 FlushDataFiles();
00537 ResetTimes();
00538 strcpy(HWConfig.fieldIdent,DEFAULT_STATION);
00539 if (manualStationName[0]!=0)
00540 {
00541 strcpy(HWConfig.fieldIdent,manualStationName);
00542 strcpy(oldHWConfig.fieldIdent,manualStationName);
00543 }
00544 strcpy(HWConfig.serialNumber,DEFAULT_SERIAL);
00545 }
00546
00547
00548 score=0;
00549 if (nTitan1+nTitan2)
00550 {
00551 score+=(100*nTitan2)/(nTitan1+nTitan2);
00552 }
00553
00554
00555 format=lastFormat;
00556 if ((score==0)&&(nTitan1>0))
00557 {
00558 format=1;
00559 }
00560 if (score>50)
00561 {
00562 format=2;
00563 }
00564
00565 if ((score<70)&&(score>5))
00566 {
00567 PrintWarning("received mixed Titan1 and Titan2 data "
00568 "frames (s=%d t1=%d, t2=%d), format=%d\n",
00569 score,nTitan1,nTitan2,format);
00570 }
00571
00572 if (format!=lastFormat)
00573 {
00574 if (lastFormat!=-1)
00575 {
00576 PrintWarning("format change detected\n");
00577
00578 CachePrintCST();
00579 FlushDataFiles();
00580 if (cutDtOnReset==0)
00581 {
00582 FlushNonDataFiles();
00583 }
00584 ResetTimes();
00585 strcpy(HWConfig.fieldIdent,DEFAULT_STATION);
00586 if (manualStationName[0]!=0)
00587 {
00588 strcpy(HWConfig.fieldIdent,manualStationName);
00589 strcpy(oldHWConfig.fieldIdent,manualStationName);
00590 }
00591 strcpy(HWConfig.serialNumber,DEFAULT_SERIAL);
00592 }
00593 lastFormat=format;
00594 if (HWConfig.formatVersion==-1)
00595 {
00596 PrintDebug("HWConfig.formatVersion=%d\n",format);
00597 HWConfig.formatVersion=format;
00598 }
00599 }
00600
00601
00602
00603 if (InWindow())
00604 {
00605
00606 extractingData=Extracting(0);
00607
00608 dataExtracted=1;
00609
00610
00611 for (i=0; i<MAX_CHANNEL; i++)
00612 {
00613
00614 if (((1LL<<i)&(auxiliaryMask|channelMask))&&(HWConfig.srObs[i].npts>=minSamplesPerFrame))
00615 {
00616
00617 if (extractTrig>0)
00618 {
00619 AddTriggerData(i);
00620 }
00621
00622 if (extractWav>0)
00623 {
00624 AddWavData(i);
00625 }
00626
00627 if (extractGif>0)
00628 {
00629 AddGifData(i);
00630 }
00631
00632 if (extractSeg2>0)
00633 {
00634 AddSeg2Data(i);
00635 }
00636
00637 if (extractSegy>0)
00638 {
00639 AddSegyData(i);
00640 }
00641
00642 if (extractSis>0)
00643 {
00644 AddSisData(i);
00645 }
00646
00647 if (extractSac>0)
00648 {
00649 AddSacData(i);
00650 }
00651
00652 if (extractMSeed>0)
00653 {
00654 AddMSeedData(i);
00655 }
00656
00657 if (extractAsc>0)
00658 {
00659 AddAscData(i);
00660 }
00661
00662 if (extractBin>0)
00663 {
00664 AddBinData(i);
00665 }
00666 }
00667
00668 if (infoHeader[i].uncorrected.usecond<=0)
00669 {
00670 memcpy(&(infoHeader[i].uncorrected),¤tTime,
00671 sizeof(struct Titan2TimeType));
00672 }
00673 }
00674 }
00675
00676
00677 memcpy(&oldHWConfig,&HWConfig,
00678 sizeof(struct HWConfigType)-sizeof(HWConfig.data));
00679
00680
00681 for (i=0; i<MAX_CHANNEL; i++)
00682 {
00683 HWConfig.srObs[i].npts=0;
00684 HWConfig.srExp[i].npts=0;
00685 lastData[i]=0;
00686 }
00687
00688 nTitan1=nTitan2=0;
00689 }
00690
00694 int64_t ProcessDataFrame(const unsigned char *frame)
00695 {
00696 int ch,i2,i5,srq,rate,n;
00697 int32_t dummy=0;
00698 int32_t d;
00699 int triplet=0;
00700 int64_t ret=0;
00701 register int i;
00702
00703
00704
00705
00706
00707
00708
00709
00710 switch (frame[11]&0x7)
00711 {
00712
00713 case DAT1_FRAME:
00714
00715 case DAT0_FRAME:
00716
00717 triplet=(frame[9]&0xf0)>>4;
00718 ret=7<<triplet;
00719
00720 if ((auxiliaryMask!=0)&&(triplet==12))
00721 {
00722 ProcessAux1Frame(frame);
00723 break;
00724 }
00725
00726 if (triplet>3)
00727 {
00728 break;
00729 }
00730
00731 nTitan1++;
00732 if (3*triplet+2>=MAX_CHANNEL)
00733 {
00734 PrintError("triplet number is %d\n",triplet);
00735 break;
00736 }
00737
00738
00739
00740
00741
00742
00743
00744 rate=frame[10]&0xf;
00745 for (i=0; i<3; i++)
00746 {
00747 int mch=0;
00748 ch=3*triplet+i;
00749 mch=ch-6;
00750 if (mch<0)
00751 {
00752 mch=ch;
00753 }
00754
00755 switch (frame[9]&0x0f)
00756 {
00757 case 8: case 9: case 10: case 11:
00758 case 12: case 13: case 14: case 15:
00759 HWConfig.srExp[ch].div=1<<(3-((frame[9]&0xf)-0x10));
00760 break;
00761 case 0:
00762 case 1:
00763 case 2:
00764 case 3:
00765 HWConfig.srExp[ch].div=1<<(3-(frame[9]&0x0f));
00766 break;
00767 default:
00768 PrintError("%d %#0x\n",
00769 frame[9]&0xf,
00770 frame[10]&0xff);
00771
00772
00773
00774
00775
00776 return(0);
00777 }
00778 if (HWConfig.srExp[ch].div<1)
00779 {
00780 HWConfig.srExp[ch].npts=0;
00781 PrintError("null divider for data frame "
00782 "(%d,%d,%d)\n",
00783 ch,mch,frame[9]&0x0f);
00784 return(0);
00785 }
00786
00787 HWConfig.srObs[ch].base=HWConfig.srExp[ch].base=HWConfig.srExp[mch].base;
00788
00789 HWConfig.srExp[ch].npts=128*HWConfig.srExp[mch].div/HWConfig.srExp[ch].div;
00790
00791 if ((rate==1)&&(HWConfig.srObs[ch].npts==0))
00792 {
00793 d=frame2i24(frame,3*i);
00794 d<<=8;
00795 d>>=8;
00796
00797 HWConfig.data[ch][HWConfig.srObs[ch].npts]=lastData[ch]+d;
00798 lastData[ch]=HWConfig.data[ch][HWConfig.srObs[ch].npts];
00799
00800 HWConfig.data[ch][HWConfig.srObs[ch].npts]+=
00801 HWConfig.absOffset[ch];
00802 if (HWConfig.srObs[ch].npts<MAX_SAMPLE-1)
00803 {
00804 HWConfig.srObs[ch].npts++;
00805 }
00806 else
00807 {
00808 PrintError("too many samples for channel"
00809 " %d, clearing buffer\n",
00810 ch);
00811 HWConfig.srObs[ch].npts=0;
00812 }
00813 }
00814
00815 else
00816 {
00817
00818 if (extracting)
00819 {
00820 dummy=frame2i24(frame,3*i);
00821 for (n=0; n<rate; n++)
00822 {
00823
00824 d=dummy&msf[rate][0];
00825
00826 dummy>>=msf[rate][1];
00827
00828 d<<=msf[rate][2];
00829 d>>=msf[rate][2];
00830 HWConfig.data[ch][HWConfig.srObs[ch].npts]=lastData[ch]+d;
00831 lastData[ch]=HWConfig.data[ch][HWConfig.srObs[ch].npts];
00832
00833 HWConfig.data[ch][HWConfig.srObs[ch].npts]+=
00834 HWConfig.absOffset[ch];
00835 HWConfig.srObs[ch].npts++;
00836 }
00837 }
00838
00839 else
00840 {
00841 HWConfig.srObs[ch].npts+=rate;
00842 }
00843 }
00844 }
00845 break;
00846 case DAT2_FRAME:
00847 if (HWConfig.formatVersion==1)
00848 {
00849 PrintDebug("DAT2 with formatVersion=1: %s\n\t%s @ %d\n",
00850 StrFrame(frame),
00851 currentFile,SafeFtell(