/****************************************************************************/
/*                                                                          */
/*                                  NOR.h                                   */
/*                                                                          */
/****************************************************************************/
/****************************************************************************/
/* This is a C language header file that defines the following routines for */
/* reading data from a NOR-format data file into memory.                    */
/* The NOR format the new data format used by Auroral Observatory of the    */
/* University of NORmsa since January 2009 in storing magnetometer data     */
/* from their stations (e.g. Svalbard stations + danish stations). The name */
/* NOR is not the official name of this data format but since I haven't     */
/* seen any name used for this data format I will call it NOR (=Norway).    */
/* These routines will handle correctly 1 second (16 bit), 10 second        */
/* and one minute data files but not 1 second (32 bit) files.               */
/*                                                                          */
/*  -------------------  Read data from a NOR-format file   --------------  */
/*                                                                          */
/* long ReadNOR(char *FileName, NetworkPtr Network, char *StartTimeStr,     */
/*               char *EndTimeStr, char *StationList, Time_sec SampleStep)  */
/*      FileName        Name of the NOR-file containing the data.           */
/*                      If FileName is nil or zero length then standard     */
/*                      input is used.                                      */
/*      Network         Structure into which the data is read.              */
/*                      Network_struct is defined in file MagnData.h.       */
/*      StartTimeStr    String defining the first time to be read from the  */
/*                      file (String delimited with '\0').                  */
/*                      String format : YYMMDDHH. If HH is missing then     */
/*                      00 is assumed. If StartTimeStr is nil or zero       */
/*                      length then start of file is assumed.               */
/*      EndTimeStr      String defining the time not to be read anymore     */
/*                      Format as before YYMMDDHH. If nil or zero           */
/*                      length then data is read until end of file.         */
/*      StationList     List of stations included in the data. The stations */
/*                      must be specified with three letter codes (e.g.SOR) */
/*                      and must be separeted with a space in the list      */
/*                      (e.g. 'SOR MAS KEV'). If nil or zero length then    */
/*                      all stations in the file are included.              */
/*      SampleStep      Time separation between successive data points.     */
/*                                                                          */
/*      The function result indicates how successfull the file reading was: */
/*          0: OK               : File read successfully                    */
/*          1: FileError        : Failed to read the given file             */
/*          2: OutOfMemory      : Memory allocation failed during reading   */
/*          3: FileFormatError  : Data file is not a NOR file               */
/*                                                                          */
/*                                                                          */
/*  ----------------------------------------------------------------------- */
/*                                                                          */
/*  NOTE:   These routines are intended to be used only in dealing with     */
/*          IMAGE magnetometer network data files in NOR-format.            */
/*          The routines defined here make assumptions about the data       */
/*          content. The routines might work also with other than IMAGE     */
/*          data files but this is not guaranteed.                          */
/*                                                                          */
/****************************************************************************/
/****************************************************************************/
/*                            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.01  12.02.2020                            */
/****************************************************************************/
/****************************************************************************/
/*  Version history:                                                        */
/*                                                                          */
/*  1.01 12.02.2020 Date strings must be YYYYMMDD (four digits for year)    */
/*                  Replaced StrToSecs -> StrToSecsC and                    */
/*                  SecsToStr -> SecsToStrC.                                */
/*  1.0  29.7.2009  First official release                                  */
/****************************************************************************/


#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <ctype.h>


/*--------------------------------------------------------------------------*/
/* Define the offsets of various fields from the start of the NOR block.    */
/* See the definition of NOR record for the explanation of various fields.  */
/*--------------------------------------------------------------------------*/

// Offsets of the various subblocks
#define NOR_Hd_block         0      // Offset of the header block
#define NOR_S1_block       117      // Offset of the scale values #1 block
#define NOR_Q1_block       249      // Offset of the QMN values #1 block
#define NOR_S2_block       315      // Offset of the scale values #2 block
#define NOR_Q2_block       447      // Offset of the QMN values #2 block
#define NOR_data_block     513      // Offset of the data block

// Offsets of the various fields within header block
#define NOR_StationNbr       0
#define NOR_StationName      1
#define NOR_StationID       22
#define NOR_Latitude        28
#define NOR_Longitude       39
#define NOR_Height          50
#define NOR_Year            52
#define NOR_Month           54
#define NOR_Day             56
#define NOR_Hour            58
#define NOR_Temp_magn       60
#define NOR_Temp_sensor     66
#define NOR_Exitation       72
#define NOR_Origin          78
#define NOR_Timing_Origin   79
#define NOR_Instrument      80
#define NOR_Filter         101
#define NOR_Timestamp      102
#define NOR_Free           103
#define NOR_Dno1           113
#define NOR_Dno2           115

// Offsets of the various fields within Scale values block
#define NOR_S_Xa             0
#define NOR_S_Ya            11
#define NOR_S_Za            22
#define NOR_S_X0            33
#define NOR_S_Y0            44
#define NOR_S_Z0            55
#define NOR_S_D0            66
#define NOR_S_I0            77
#define NOR_S_N1            88
#define NOR_S_N2            99
#define NOR_S_N3           110
#define NOR_S_N4           121

// Offsets of the various fields within QMNvalues block
#define NOR_Q_QmX            0
#define NOR_Q_QmY           11
#define NOR_Q_QmZ           22
#define NOR_Q_QmD           33
#define NOR_Q_QmI           44
#define NOR_Q_QmF           55

#define MissingNOR      -32768
#define Deg_to_Rad      0.0174532925

typedef char *NORptr;


/*------------------------------------------------------*/
/* The data structure for single NOR binary data block  */
/*------------------------------------------------------*/

struct NORdatablock_struct {
    Time_sec    StartTime;      /* Time of the first data point */
    char        *Data;          /* Pointer to data area         */
    long        status;
};

typedef struct NORdatablock NORdatablock_struct;
typedef NORdatablock_struct *NORdatablockPtr;


/*--------------------------------------------------------------------------*/
/* Function prototypes                                                      */
/*--------------------------------------------------------------------------*/

static void P_to_C_str(char *P_str, char *C_str);
static double atof_P_str(char *P_str);
static long GetWordLH(char *p);
static long GetShortLH(char *p);
static Time_sec GetNORTime(NORptr NOR);
static long ReadNORBlock(FILE *NORfile, NORptr NOR, long BlockSize);
static void GetNORStationID(char *StationID,NORptr NOR);
static void SetStationParamsNOR(StationPtr Station,NORptr NOR, Time_sec StartTime, long SampleStep);
static void FindScaleValues(NORptr Block, double *Xscale, double *Yscale, double *Zscale,
                            double *D0, double *I0, double *cosI0, double *sinI0);
static void GetNORData(StationPtr Station,NORptr NOR);
long ReadNOR(char *FileName, NetworkPtr Network, char *StartTimeStr,
              char *EndTimeStr, char *StationList, Time_sec SampleStep);


/*==========================================================================*/
/*              Routines for reading NOR-format files                       */
/*==========================================================================*/

/*--------------------------------------------------------------------------*/
/*  Convert a Pascal string into a C string                                 */
/*--------------------------------------------------------------------------*/

static void P_to_C_str(char *P_str, char *C_str)
{
    register long i,count;
    register char *p = P_str+1;
    register char *q = C_str;

    count = *P_str;
    if (count > 10) count = 10;     // Maximum length = 10
    for (i=0;i<count;i++) *q++ = *p++;
    *q = '\0';
}

/*--------------------------------------------------------------------------*/
/*  Convert a Pascal string representing a floating point number into double*/
/*--------------------------------------------------------------------------*/

static double atof_P_str(char *P_str)
{
    char C_str[20];

    P_to_C_str(P_str,C_str);
    return (atof(C_str));
}


/*--------------------------------------------------------------------------*/
/*  Get a two-byte binary number. The order of the bytes is : Low-High      */
/*  We cannot read the binaries with (short) p, where p is defined as       */
/*  (short) *p, since some processors use high-low byte ordering. In the    */
/*  NOR format the ordering is defined as low-high. We could check the      */
/*  type of the processor and then swap the bytes if necessary but that     */
/*  wouldn't be any simpler.                                                */
/*--------------------------------------------------------------------------*/

static long GetWordLH(char *p)      // Get unsigned integer
{
    register long high,low;

    low  = 0x000000FF & (*p++);
    high = 0x000000FF & (*p);
    return((high << 8) | low);
}

static long GetShortLH(char *p)     // Get signed integer
{
    register long high,low;

    low  = 0x000000FF & (*p++);
    high = (signed char) (*p);
    return((high << 8) | low);
}


/*--------------------------------------------------------------------------*/
/*  Get the start time in seconds of a NOR data block.                      */
/*--------------------------------------------------------------------------*/

static Time_sec GetNORTime(NORptr NOR)
{
    Time_struct MyTime;

    MyTime.tm_year  = GetWordLH(NOR+NOR_Year)-1900;
    MyTime.tm_mon   = GetWordLH(NOR+NOR_Month)-1;
    MyTime.tm_mday  = GetWordLH(NOR+NOR_Day);
    MyTime.tm_hour  = GetWordLH(NOR+NOR_Hour);
    MyTime.tm_min   = 0;
    MyTime.tm_sec   = 0;

// fprintf(stderr,"%04d %02d %02d   %02d\n",MyTime.tm_year+1900,MyTime.tm_mon+1,MyTime.tm_mday,MyTime.tm_hour);

    return TmToSecs(&MyTime);
}


/*--------------------------------------------------------------------------*/
/* Read one data block from the NOR-format file into NOR block. Function    */
/* returns 1 if reading was successful and 0 if the end of file is reached. */
/*--------------------------------------------------------------------------*/

static long ReadNORBlock(FILE *NORfile, NORptr NOR, long BlockSize)
{
    long FileStatus;

    FileStatus = fread(NOR,BlockSize,1,NORfile);
    if (FileStatus == 0)
        return (0);
    else
        return (1);
}


/*--------------------------------------------------------------------------*/
/*  Find out the three letter station ID from the ID-number code in the     */
/*  NOR block.                                                              */
/*--------------------------------------------------------------------------*/

static void GetNORStationID(char *StationID,NORptr NOR)
{
    char ID[20];
    int i;

//  Before 2011 the station ID was determined from one byte station number.
//  With inclusion of Danish magnetometer data the system changed and the
//  station number was put to 99 (= missing station) for all stations.
//     ID = *(NOR+NOR_Hd_block+NOR_StationNbr);
//     while (i<IMAGEStationCount) {
//         if (StationInfoTable[i++].TRO_ID == ID) break;
//     }
//     strcpy(StationID,StationInfoTable[i-1].StationID);

    P_to_C_str(NOR+NOR_Hd_block+NOR_StationID,ID);
    for (i=0;i<3;i++) StationID[i] = toupper(ID[i]);
    StationID[3] = '\0';
}


/*--------------------------------------------------------------------------*/
/*  Fill some fields in the Station_struct by copying them from the NOR     */
/*  block. AddNewStation procedure in the MagnData.h file will set other    */
/*  station parameters.                                                     */
/*--------------------------------------------------------------------------*/

static void SetStationParamsNOR(StationPtr Station,NORptr NOR, Time_sec StartTime, long SampleStep)
{
    char Lat_str[10],Lon_str[10];

    /* --- Set the three letter station ID --- */
    GetNORStationID(Station->StationID,NOR);
    Station->StationID[3] = '\0';

    /* --- Set the coordinates of the station, unit = 0.01 degrees --- */
    P_to_C_str(NOR+NOR_Latitude,Lat_str);
    P_to_C_str(NOR+NOR_Longitude,Lon_str);

    Station->Latitude = RoundFloat(100*atof(Lat_str));
    Station->Longitude = RoundFloat(100*atof(Lon_str));

    /* Set the component markers */
    strcpy(Station->Components,"XYZ");

    /* --- Set StartTime and TimeStep, GetNORData-routine sets EndTime --- */
    Station->StartTime = StartTime;
    Station->EndTime   = Station->StartTime;
    Station->TimeStep  = SampleStep;
}


/*--------------------------------------------------------------------------*/
/*  Find the scale values , offsets and declination and inclination values  */
/*  from the Scale values data block.                                       */
/*--------------------------------------------------------------------------*/

static void FindScaleValues(NORptr Block, double *Xscale, double *Yscale, double *Zscale,
                            double *D0, double *I0, double *cosI0, double *sinI0)
{
    Xscale[0] = atof_P_str(Block + NOR_S_Xa);   // index 0: scale values
    Yscale[0] = atof_P_str(Block + NOR_S_Ya);
    Zscale[0] = atof_P_str(Block + NOR_S_Za);
    Xscale[1] = atof_P_str(Block + NOR_S_X0);   // index 1: offsets
    Yscale[1] = atof_P_str(Block + NOR_S_Y0);
    Zscale[1] = atof_P_str(Block + NOR_S_Z0);
    *D0 = atof_P_str(Block + NOR_S_D0);
    *I0 = atof_P_str(Block + NOR_S_I0);
    *sinI0 = sin((*I0)*Deg_to_Rad);             // D0 and I0 are in degrees
    *cosI0 = cos((*I0)*Deg_to_Rad);
}


/*--------------------------------------------------------------------------*/
/*  Copy data from NOR block into proper OneDayBlock in Station structure.  */
/*--------------------------------------------------------------------------*/

static void GetNORData(StationPtr Station,NORptr NOR)
{
    OneDayPtr p;
    long *Xptr, *Yptr, *Zptr, *Tptr;
    double Xscale[2];   // scales and offsets. 0 = scale, 1 = offset
    double Yscale[2];
    double Zscale[2];
    double D0, I0, cosI0, sinI0;
    double XX, YY, ZZ, tx, H, D;
    long Scale_2_start;
    long Count;
    long DataOrigin;
    long i, step;
    long X,Y,Z;

    /* Find the last OneDayBlock */
    for (p = Station->FirstDay; p->NextDay != NULL; p = p->NextDay);

    /* Set the temperature value */
    Tptr  = (p->T)+(p->TCount);
    *Tptr = GetShortLH(NOR + NOR_Temp_magn);
    if (*Tptr == MissingNOR) *Tptr = MissingValue;
    p->TCount++;

    /* Get the baselines and scalefactors for the first block */
    FindScaleValues(NOR+NOR_S1_block, Xscale, Yscale, Zscale, &D0, &I0, &cosI0, &sinI0);

    Scale_2_start = GetWordLH(NOR + NOR_Dno1);  // Start index of second scale value set


    /* Get the number of datapoints */
    switch (Station->TimeStep) {
        case 1  : Count = 10800; break; // 1 second short
        case 10 : Count =  1080; break; // 10 seconds
        case 60 : Count =  1440; break; // 1 minute
        default : Count =  1080;        // 10 seconds
    }

    /* Get raw data, multiply by scale coefficients, add offsets */
    /* and rotate into geographic coordinates if necessary.      */

    DataOrigin = *(NOR + NOR_Origin);   // 0 = raw data, 1 = geographic XYZ

    Xptr = (p->X)+(p->XCount);
    Yptr = (p->Y)+(p->YCount);
    Zptr = (p->Z)+(p->ZCount);
    step = 0;

    for (i=0;i<Count;i++) {
        X = GetShortLH(NOR + NOR_data_block + step);
        Y = GetShortLH(NOR + NOR_data_block + 2*Count + step);
        Z = GetShortLH(NOR + NOR_data_block + 4*Count + step);

        if ((X == MissingNOR) || (Y == MissingNOR) || (Z == MissingNOR)) {
            *Xptr++ = MissingValue;
            *Yptr++ = MissingValue;
            *Zptr++ = MissingValue;
        }
        else {
            XX = Xscale[0]*X + Xscale[1];
            YY = Yscale[0]*Y + Yscale[1];
            ZZ = Zscale[0]*Z + Zscale[1];

            if (DataOrigin == 0) {      // Raw data, rotate to geographic coordinates
                tx = ZZ*cosI0 + XX*sinI0;           // temporary 'X'
                H  = sqrt(tx*tx + YY*YY);
                D  = atan2(YY,tx) + Deg_to_Rad*D0;  // Declination in radians
                *Xptr++ = RoundFloat(10.0*H*cos(D));                // X geographic
                *Yptr++ = RoundFloat(10.0*H*sin(D));                // Y geographic
                *Zptr++ = RoundFloat(10.0*(ZZ*sinI0 - XX*cosI0));   // Z geographic
            }
            else {                      // already geographic XYZ
                *Xptr++ = RoundFloat(10.0*XX);      // Unit = 0.1 nT
                *Yptr++ = RoundFloat(10.0*YY);
                *Zptr++ = RoundFloat(10.0*ZZ);
            }
        }

        if (i == (Scale_2_start-1)) {   // Check if we need to update the set of coefficients
            FindScaleValues(NOR+NOR_S2_block, Xscale, Yscale, Zscale, &D0, &I0, &cosI0, &sinI0);
        }
        step += 2;
    }

    /* --- Update counters --- */
    p->XCount += Count;
    p->YCount += Count;
    p->ZCount += Count;

    /* --- Set the end time --- */
    Station->EndTime = GetNORTime(NOR)+Count*Station->TimeStep;
}


/*--------------------------------------------------------------------------*/
/*  Copy data from NOR block into proper OneDayBlock in Station structure.  */
/*--------------------------------------------------------------------------*/

static void ReplaceNORData(StationPtr Station, Time_sec CurrTime, NORptr NOR)
{
    OneDayPtr p;
    long *Xptr, *Yptr, *Zptr, *Tptr;
    double Xscale[2];   // scales and offsets. 0 = scale, 1 = offset
    double Yscale[2];
    double Zscale[2];
    double D0, I0, cosI0, sinI0;
    double XX, YY, ZZ, tx, H, D;
    long Scale_2_start;
    long Count;
    long DataOrigin;
    long i, step;
    long X,Y,Z;
    long TCount, DataCount;
    Time_sec t;

    /* Find the correct OneDayBlock */
    p = Station->FirstDay;
    t = Station->StartTime;
    while (CurrTime >= t + 86400) {
        t += 86400;
        p = p->NextDay;
    }

    TCount = (CurrTime - t)/(3*3600);   // One temp value for each three hours
    /* Set the temperature value */
    Tptr  = (p->T)+(TCount);
    *Tptr = GetShortLH(NOR + NOR_Temp_magn);
    if (*Tptr == MissingNOR) *Tptr = MissingValue;

    /* Get the baselines and scalefactors for the first block */
    FindScaleValues(NOR+NOR_S1_block, Xscale, Yscale, Zscale, &D0, &I0, &cosI0, &sinI0);

    Scale_2_start = GetWordLH(NOR + NOR_Dno1);  // Start index of second scale value set

    /* Get the number of datapoints */
    switch (Station->TimeStep) {
        case 1  : Count = 10800; break; // 1 second short
        case 10 : Count =  1080; break; // 10 seconds
        case 60 : Count =  1440; break; // 1 minute
        default : Count =  1080;        // 10 seconds
    }

    /* Get raw data, multiply by scale coefficients, add offsets */
    /* and rotate into geographic coordinates if necessary.      */

    DataOrigin = *(NOR + NOR_Origin);   // 0 = raw data, 1 = geographic XYZ

    DataCount = Count*(CurrTime - t)/(3*3600);  // Three hour blocks

// fprintf(stderr,"%d %d\n", TCount, DataCount);

    Xptr = (p->X)+(DataCount);
    Yptr = (p->Y)+(DataCount);
    Zptr = (p->Z)+(DataCount);
    step = 0;

    for (i=0;i<Count;i++) {
        X = GetShortLH(NOR + NOR_data_block + step);
        Y = GetShortLH(NOR + NOR_data_block + 2*Count + step);
        Z = GetShortLH(NOR + NOR_data_block + 4*Count + step);

        if ((X == MissingNOR) || (Y == MissingNOR) || (Z == MissingNOR)) {
            *Xptr++ = MissingValue;
            *Yptr++ = MissingValue;
            *Zptr++ = MissingValue;
        }
        else {
            XX = Xscale[0]*X + Xscale[1];
            YY = Yscale[0]*Y + Yscale[1];
            ZZ = Zscale[0]*Z + Zscale[1];

            if (DataOrigin == 0) {      // Raw data, rotate to geographic coordinates
                tx = ZZ*cosI0 + XX*sinI0;           // temporary 'X'
                H  = sqrt(tx*tx + YY*YY);
                D  = atan2(YY,tx) + Deg_to_Rad*D0;  // Declination in radians
                *Xptr++ = RoundFloat(10.0*H*cos(D));                // X geographic
                *Yptr++ = RoundFloat(10.0*H*sin(D));                // Y geographic
                *Zptr++ = RoundFloat(10.0*(ZZ*sinI0 - XX*cosI0));   // Z geographic
            }
            else {                      // already geographic XYZ
                *Xptr++ = RoundFloat(10.0*XX);      // Unit = 0.1 nT
                *Yptr++ = RoundFloat(10.0*YY);
                *Zptr++ = RoundFloat(10.0*ZZ);
            }
        }

        if (i == (Scale_2_start-1)) {   // Check if we need to update the set of coefficients
            FindScaleValues(NOR+NOR_S2_block, Xscale, Yscale, Zscale, &D0, &I0, &cosI0, &sinI0);
        }
        step += 2;
    }
}



/*--------------------------------------------------------------------------*/
/* Add missing data block                                                   */
/*--------------------------------------------------------------------------*/

static void AddMissingNORBlock(StationPtr Station, long Count)
{
    OneDay_struct *p;

    /* Find the last OneDayBlock */
    for (p = Station->FirstDay; p->NextDay != NULL; p = p->NextDay);

    /* We only have to update counters since the data areas are already */
    /* initialized with missing data                                    */
    p->XCount += Count;
    p->YCount += Count;
    p->ZCount += Count;
    p->TCount++;
}



/*--------------------------------------------------------------------------*/
/*  Here is the routine for reading NOR-format data file into memory.       */
/*  The data is read into a Network structure (see MagnData.h for details). */
/*  See the comments at the beginning of this file for more details about   */
/*  the parameters for this routine.                                        */
/*                                                                          */
/*  The rounite works by reading the data file one block at a time and adds */
/*  that data in the particular stations data block.                        */
/*  For each new station that is encountered space for one day is allocated.*/
/*  If the day becomes full a new one day block is allocated. At the end of */
/*  the routine all one day blocks are combined into a single data block.   */
/*--------------------------------------------------------------------------*/

long ReadNOR(char *FileName, NetworkPtr Network, char *StartTimeStr,
              char *EndTimeStr, char *StationList, Time_sec SampleStep)
{
    FILE    *NORfile;               /* The data file to be processed        */
    char    *NORblock;              /* Pointer to data record               */
    Time_sec CurrTime;              /* Time of the current data block       */
    Time_sec PrevTime = MIN_TIME;   /* Time of the previous data block      */
    Time_sec StartTime;             /* Time of first record to be read      */
    Time_sec EndTime;               /* Time of record not read anymore      */
    StationPtr Station;             /* Pointer to the current station       */
    char    StationID[4];           /* Station ID of the current station    */
    long    Count;                  /* Number of data points in one block   */
    long    BlockStatus;            /* Status number for block reading      */
    long    BlockSize;              /* Size of a single data block          */
    InitNetwork(Network);           /* Initialize the fields of Network     */

    struct NORdatablock_struct Blocks[300];     /* One block = 3 hours, One month = 31*8 = 248 blocks */
    long BlockCount = 0;
    long BlockNo = 0;
    long i;



    char TimeStr[20];

    /* --- Set the start and end times --- */
    if ((StartTimeStr == NULL) || (StartTimeStr[0] == '\0'))
        StartTime = MIN_TIME;
    else StartTime = StrToSecsC(StartTimeStr,0);

    if ((EndTimeStr == NULL) || (EndTimeStr[0] == '\0'))
        EndTime = MAX_TIME-(Time_sec)86400;
    else EndTime = StrToSecsC(EndTimeStr,0);

    /* --- Allocate space for data record --- */
    switch (SampleStep) {
        case 1  : BlockSize = NOR_data_block + 64800; Count = 10800; break; // 1 second short
        case 10 : BlockSize = NOR_data_block + 6480;  Count =  1080; break; // 10 seconds
        case 60 : BlockSize = NOR_data_block + 8640;  Count =  1440; break; // 1 minute
        default : BlockSize = NOR_data_block + 6480;  Count =  1080;        // 10 seconds
    }

    NORblock = (char *) malloc(BlockSize);


    /* --- Try to open the file --- */
    if ((FileName == NULL) || (FileName[0] == '\0')) NORfile = stdin;
    else if ((NORfile = fopen(FileName, "rb")) == NULL) return(FileError);


    // Read all data blocks into memory. If same time appears twice then replace the first
    // data with the second one.

     while (BlockStatus = ReadNORBlock(NORfile,NORblock,BlockSize)) {
        // Check the three letter stationID. If this block is garbage (happens sometimes) then
        // stationID is garbage also.
        GetNORStationID(StationID,NORblock);
        if ((Network->StationList != NULL) && (strncmp(StationID, Network->StationList->StationID,3) != 0)) {
            continue;   // This block is for other station or garbage, so neglect it
        }

        CurrTime=GetNORTime(NORblock);

        BlockNo = BlockCount;
        for (i=0; i < BlockCount; i++) {
            if (Blocks[i].StartTime == CurrTime) {
                BlockNo = i;
                break;
            }
        }

        Blocks[BlockNo].StartTime = CurrTime;
        Blocks[BlockNo].Data = (char *) malloc(BlockSize);
        memcpy(Blocks[BlockNo].Data, NORblock, BlockSize);
        Blocks[BlockNo].status = 0;

        if (BlockNo == BlockCount) BlockCount++;
     }



    /* --- Read data into one day blocks ---*/
    while (1) {

        // Get the block with earliest time which has not been read already

        CurrTime = MAX_TIME;
        BlockNo = -1;
        for (i=0; i < BlockCount; i++) {
            if ((Blocks[i].status == 0) && (Blocks[i].StartTime <= CurrTime)) {
                CurrTime = Blocks[i].StartTime;
                BlockNo = i;
            }
        }

        if (BlockNo == -1) break;

// fprintf(stderr,"%ld\n", CurrTime);

        Blocks[BlockNo].status = 1;

        if (CurrTime >= EndTime+(Time_sec)86400) break;

        GetNORStationID(StationID, Blocks[BlockNo].Data);

        if ((CurrTime >= PrevTime) && (CurrTime >= StartTime) && (CurrTime < EndTime)
            && StationInList(StationID,StationList))
        {
            Station = FindStation(Network,StationID);
            if (Station == NULL) {      /* First occurrence of this station */
                 Station = AddNewStation(Network,SampleStep,Count*SampleStep);
                if (Station != NULL) {
                    SetStationParamsNOR(Station,Blocks[BlockNo].Data,CurrTime,SampleStep);
                }
                else        /* Unable to allocate memory for station data   */
                    return(OutOfMemory);
            }

            /* --- Add possible missing data blocks ---                */
            /* PrevTime is the expected time of the current data block.*/
            /* If it is less than the time read from the block then    */
            /* add missing data blocks                                 */

            PrevTime = Station->EndTime;
            while (PrevTime < CurrTime) {
                if (OneDayFull(Station,'Z')) {
                    if (AddOneDay(Station))
                        return(OutOfMemory);
                }

                AddMissingNORBlock(Station, Count);
                PrevTime += Count*SampleStep;
            }

            if (OneDayFull(Station,'Z')) {
                if (AddOneDay(Station))
                    return(OutOfMemory);    /* Unable to allocate memory    */
            }
            GetNORData(Station, Blocks[BlockNo].Data);
        }
    }
    fclose(NORfile);

    /* --- Combine one day blocks into one large block ---*/
    Station = Network->StationList;
    while (Station != NULL) {
        if (CombineData(Station) == OutOfMemory) return(OutOfMemory);
        Station = Station->Next;
    }

    return OK;
}