/****************************************************************************/
/*                                                                          */
/*                              Find_info.c                                 */
/*                                                                          */
/****************************************************************************/
/****************************************************************************/
/* This program reads IAGA, GADF, Dump or WDC format data files containing  */
/* magnetometer data and writes information about data content (Station     */
/* names, start and end times, data availability percent) into standard     */
/* output                                                                   */
/*                                                                          */
/* Usage:                                                                   */
/*    Find_info [-f] [-p] [<] File                                          */
/*      [-f Dataformat]     Format of the data file. Possible values for    */
/*                          Dataformat are IAGA (default),GADF,Dump or WDC. */
/*      [-p]                Prevents usage-text from appearing. Useful when */
/*                          reading data from standard input.               */
/*      File                Name of data file.                              */
/*                                                                          */
/****************************************************************************/
/****************************************************************************/
/*                            Lasse Hakkinen                                */
/*                    Finnish Meteorological Institute                      */
/*                        Department of Geophysics                          */
/*                              P.O.Box 503                                 */
/*                      FIN-00101, Helsinki, Finland                        */
/*                      e-mail: Lasse.Hakkinen@fmi.fi                       */
/*                      phone : (+358)-9-19294634                           */
/*                      fax   : (+358)-9-19294603                           */
/*                                                                          */
/*                      version 1.03        22.05.2001                      */
/****************************************************************************/
/****************************************************************************/
/*  Version history:                                                        */
/*                                                                          */
/*  1.03 22.05.2001 Changed the output format by adding the time of Min(X)  */
/*                  and Max(|dB/dt|) and time of Max(|dX/dt|).              */
/*  1.02 02.07.1999 Fixed a Y2K bug in NewTime.h file which resulted in     */
/*                  incorrect year in date strings if year >= 2000.         */
/*  1.01 21.10.1998 Accommodated to the changes in Dump.h (0.1 nT to        */
/*                  0.01 nT resolution).                                    */
/*  1.0  04.06.1998 First official release                                  */
/****************************************************************************/

#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include "Usage.h"
#include "MagnData.h"
#include "IAGA.h"
#include "GADF.h"
#include "WDC.h"
#include "Dump.h"

char *version = "1.03";
char *date = "22.05.2001";

#define HeaderLineCount 5
#define UsageLineCount  6

char *HeaderText[HeaderLineCount] = {
"**************************************************************************",
"This program reads a magnetometer data file and writes out information    ",
"about the data content of the file (station names, start and end times,   ",
"data availability percent etc.).                                          ",
"**************************************************************************",
};

char *UsageText[UsageLineCount] = {
" [-f] [-p] [<] File                                                        ",
"      [-f Dataformat]      Format of the data file. Possible values for    ",
"                           Dataformat are IAGA (default),GADF,Dump or WDC. ",
"      [-p]                 Prevents usage text from appearing. Useful      ",
"                           when reading from standard input.               ",
"      File                 Name of data file.                              ",
};



/*--------------------------------------------------------------------------*/
/*  Find the number of missing data points in one component                 */
/*--------------------------------------------------------------------------*/

long FindMissingCount(StationPtr Station,char Comp)
{
    register long i;
    long MissCount = 0;
    long TotalCount;
    long *DataPtr;
    
    switch (Comp) {
        case 'X' : DataPtr = Station->X; break;
        case 'Y' : DataPtr = Station->Y; break;
        case 'Z' : DataPtr = Station->Z; break;
    }
    
    TotalCount = (Station->EndTime - Station->StartTime)/(Station->TimeStep);
    for (i=0;i<TotalCount;i++) {
        if (*DataPtr == MissingValue) MissCount++;
        DataPtr++;
    }
    return (MissCount);
}


/*--------------------------------------------------------------------------*/
/*  Find the times of maximum and minimum field values over speficied time  */
/*  period station, component and time.                                     */
/*--------------------------------------------------------------------------*/

void FindMaxMinTime(StationPtr Station,char Comp,Time_sec Time,long Length,
                Time_sec *MaxTime, Time_sec *MinTime)
{
    long i,Count,Value;
    long Good = 0;
    long Max,Min;
    Time_sec Tmin,Tmax;
    
    if (Length == 0)
        Count = 1;  /* Max and Min are same */
    else
        Count = Length/Station->TimeStep;

    for (i=0;i<Count;i++) {
        Value = GetDataValue(Station,Time,Comp);
        if (Value != MissingValue) {
            if (Good == 0) {    /* Is this the first good data point ? */
                Max = Value;
                Min = Value;
                Tmax = Time;
                Tmin = Time;
                Good = 1;
            }
            
            if (Value > Max) {
                Max = Value;
                Tmax = Time;
            }
            else if (Value < Min) {
                Min = Value;
                Tmin = Time;
            }
        }
        Time += Station->TimeStep;
    }
    if (Good == 0) {
        *MaxTime = MIN_TIME;
        *MinTime = MIN_TIME;
    } else {
        *MaxTime = Tmax;
        *MinTime = Tmin;
    }
}

/*--------------------------------------------------------------------------*/
/*  Find the maximum of the derivative of the field and corresponding time  */
/*  over speficied time period station, component and time.                 */
/*--------------------------------------------------------------------------*/

void FindMaxDer(StationPtr Station,char Comp,Time_sec Time,long Length,
                double *MaxValue, Time_sec *MaxTime)
{
    long i,Count,Value1,Value2;
    long Good = 0;
    double Max;
    double Der;
    Time_sec Tmax;
    
    if (Length == 0)
        Count = 1;  /* Max and Min are same */
    else
        Count = Length/Station->TimeStep;

    for (i=0;i<Count-1;i++) {
        Value1 = GetDataValue(Station,Time,Comp);
        Value2 = GetDataValue(Station,Time+Station->TimeStep,Comp);
        if ((Value1 != MissingValue) && (Value2 != MissingValue)) {
            if (Good == 0) {    /* Is this the first good data pair ? */
                Max = ((double) abs(Value2-Value1))/Station->TimeStep;
                Tmax = Time;
                Good = 1;
            }
            
            Der = ((double) abs(Value2-Value1))/Station->TimeStep;
            if (Der > Max) {
                Max = Der;
                Tmax = Time;
            }
        }
        Time += Station->TimeStep;
    }
    if (Good == 0) {
        *MaxTime = MIN_TIME;
        *MaxValue = MissingValue;
    } else {
        *MaxTime = Tmax;
        *MaxValue = Max/10.0;
    }
}


/*--------------------------------------------------------------------------*/
/*                          The main procedure                              */
/*--------------------------------------------------------------------------*/

int main(int argc, char *argv[])
{
    long    params,i;
    long    status = 0;
    long    FileCount = 0;
    char    FileName[100] = "";
    Network_struct NETWORK;             /* Structure containing the data    */
                                        /* for all stations.                */
    StationPtr  Station;                /* Pointer to the current station   */
    char FormatStr[10]      = "IAGA";   /* Format of magnetometer data file */
    char StartTimeStr[20];
    char EndTimeStr[20];
    char CenturyStr[3];
    char Comp;
    long StationNbr = 0;
    long MissingCount;
    long TotalCount;
    long Max,Min;
    Time_sec MaxTime,MinTime;
    double MaxDer;


    /*==========================*/
    /* Analyse the command line */
    /*==========================*/
    if (argc == 1) {        /* No arguments, show the usage text */
        PrintUsage(argv[0],0,HeaderLineCount,UsageLineCount);
        return OK;
    }

    for (params = 1; params < argc; params++) {
        if (*argv[params] != '-') {
            strcpy(FileName,argv[params]);
            FileCount++;
        } else {
            switch (*(argv[params]+1)) {
                case 'f' : strcpy(FormatStr,argv[++params]); break;
                case 'p' : /* do nothing */                  break;
                default  :
                    fprintf(stderr,"\n### %s: \"%s\" is not an option.\n\n",
                            argv[0], argv[params]);
                    PrintUsage(argv[0],1,HeaderLineCount,UsageLineCount);
                    return FAIL;
                    break;
            }
        }
    }


    /*=================================================*/
    /* Check that all necessary parameters are defined */
    /*=================================================*/
    
    if (FileCount > 1) {
        PrintUsage(argv[0],1,HeaderLineCount,UsageLineCount);
        return FAIL;
    }


    /*======================================*/
    /* Read data from data file into memory */
    /*======================================*/

    if (!strcmp(FormatStr,"IAGA")) {
        status = ReadIAGA(FileName,&NETWORK,NULL,NULL,NULL);
    } else
    if (!strcmp(FormatStr,"GADF")) {
        status = ReadGADF(FileName,&NETWORK,NULL,NULL,NULL);
    } else
    if (!strcmp(FormatStr,"Dump")) {
        status = ReadDump(FileName,&NETWORK,NULL,NULL,NULL);
    } else
    if (!strcmp(FormatStr,"WDC")) {
        status = ReadWDC(FileName,&NETWORK,NULL,NULL,NULL);
    } else {
        fprintf(stderr,"Illegal format type: %s\n",FormatStr);
        return FAIL;
    }

    if (status != 0) {
        if (status == FileError)
            fprintf(stderr,"### Error in opening data file ");
        if (status == OutOfMemory)
            fprintf(stderr,"### Out of memory while reading data file");
        if (status == FileFormatError)
            fprintf(stderr,"### Wrong file format in input file ");
        fprintf(stderr,"%s\n",FileName);
        return FAIL;
    }

    if (!strcmp(FormatStr,"Dump"))
        ChangeUnitsDOWN(&NETWORK);


    /*======================================*/
    /* Write header                         */
    /*======================================*/

    Station = NETWORK.StationList;
    SecsToStr(Station->StartTime,StartTimeStr);
    SecsToStr(Station->EndTime,EndTimeStr);
    
    if (Station->StartTime < YEAR2000) strcpy(CenturyStr,"19"); else strcpy(CenturyStr,"20");
    printf("\nStart time: %.2s",CenturyStr);
    printf("%.2s-%.2s-%.2s ",StartTimeStr  ,StartTimeStr+2,StartTimeStr+4);
    printf("%.2s:%.2s:%.2s ",StartTimeStr+6,StartTimeStr+8,StartTimeStr+10);

    if (Station->EndTime < YEAR2000) strcpy(CenturyStr,"19"); else strcpy(CenturyStr,"20");
    printf("\nEnd time  : %.2s",CenturyStr);
    printf("%.2s-%.2s-%.2s ",EndTimeStr  ,EndTimeStr+2,EndTimeStr+4);
    printf("%.2s:%.2s:%.2s ",EndTimeStr+6,EndTimeStr+8,EndTimeStr+10);

    printf("\nTime resolution: %d seconds",Station->TimeStep);

    printf("\n\n");
    
    printf("-------------------------------------------------------------------------------\n");
    printf("No. ID  Data avail     Max - Min (nT)    Time of   Max |dB/dt| (nT/s)  Time of\n");
    printf("           (%%)        X      Y      Z     Min X      X     Y     Z   Max|dX/dt|\n");
    printf("-------------------------------------------------------------------------------\n");


    /*======================================*/
    /* Go through all stations in NETWORK   */
    /*======================================*/

    while(Station != NULL) {

        StationNbr++;

        /*** Write station ID ***/
        
        printf("%2d: %s  ",StationNbr,Station->StationID);
        

        /*** Find the number of missing data points in each component ***/
        
        MissingCount = FindMissingCount(Station,'X');
        MissingCount += FindMissingCount(Station,'Y');
        MissingCount += FindMissingCount(Station,'Z');
        TotalCount = 3*(Station->EndTime - Station->StartTime)/(Station->TimeStep);
        printf("%6.2f  ",100.0-(100.0*MissingCount)/TotalCount);


        /*** Find the amplitudes (= Max-Min) ***/
        
        for (Comp='X';Comp <= 'Z';Comp++) {
            FindMaxMin(Station,Comp,Station->StartTime,
                Station->EndTime - Station->StartTime,&Max,&Min);
            if ((Max-Min) > 0)
                printf(" %6.1f",(Max-Min)/10.0);
            else
                printf("    ---");
        }


        /*** Find the time of Min(X) ***/

        FindMaxMinTime(Station,'X',Station->StartTime,Station->EndTime - Station->StartTime,&MaxTime,&MinTime);
        SecsToStr(MinTime,StartTimeStr);
        if (MinTime != MIN_TIME)
            printf("  %.2s:%.2s:%.2s ",StartTimeStr+6,StartTimeStr+8,StartTimeStr+10);
        else
            printf("    ---    ");


        /*** Find the maxima of time derivatives of each component ***/

        for (Comp='X';Comp <= 'Z';Comp++) {
            FindMaxDer(Station,Comp,Station->StartTime,Station->EndTime - Station->StartTime,&MaxDer,&MaxTime);
            if (MaxDer != MissingValue)
                printf(" %5.1f",MaxDer);
            else
                printf("   ---");
            if (Comp == 'X') MinTime = MaxTime;
        }

        /*** Print the time of Max |dX/dt| ***/
        
        SecsToStr(MinTime,StartTimeStr);
        if (MinTime != MIN_TIME)
            printf("   %.2s:%.2s:%.2s",StartTimeStr+6,StartTimeStr+8,StartTimeStr+10);
        else
            printf("     ---   ");


        printf("\n");


        Station = Station->Next;
    }                               /* End of Station Loop */

    printf("-------------------------------------------------------------------------------\n\n");

    FreeNetwork(&NETWORK);
    return OK;
}