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