(file) Return to jretile_manytofew.c CVS log (file) (dir) Up to [Development] / JSOC / proj / globalhs / apps

File: [Development] / JSOC / proj / globalhs / apps / jretile_manytofew.c (download)
Revision: 1.3, Tue Sep 18 15:03:27 2012 UTC (5 years, 10 months ago) by tplarson
Branch: MAIN
Changes since 1.2: +28 -5 lines
added checking of CALVER64

#define MODES(L)        ((((L)+1)*(L))/2)
#define MINIMUM(X,Y)    (((X)<(Y))?(X):(Y))
#define MAXIMUM(X,Y)    (((X)>(Y))?(X):(Y))

int jretile_manytofew(void)
{
  int newstat = 0;
  int status = DRMS_SUCCESS;
  int fetchstat = DRMS_SUCCESS;
  DRMS_RecChunking_t chunkstat = kRecChunking_None;

  char *inrecquery = NULL;
  char *outseries = NULL;
  char *segnamein = NULL;
  char *segnameout = NULL;
  DRMS_RecordSet_t *inrecset = NULL;
  DRMS_RecordSet_t *outrecset = NULL;
  int irecin, irecout, nrecsin=0, nrecsout=0, nlchunks;
  DRMS_Record_t *inrec = NULL;
  DRMS_Record_t *outrec = NULL;
  DRMS_Segment_t *segin = NULL;
  DRMS_Segment_t *segout = NULL;
  DRMS_Array_t *inarr = NULL;
  DRMS_Array_t *outarr = NULL;
  DRMS_RecLifetime_t lifetime;
  DRMS_Type_t usetype = DRMS_TYPE_FLOAT;
  int length[2], startind[2], endind[2], totallength[2];
  float *inptr, *outptr;
  long long histrecnum=-1;
  int quality;
  int mapmmax=-1;
  int sinbdivs=-1;
  double cadence=0;

  TIME tnow, UNIX_epoch = -220924792.000; /* 1970.01.01_00:00:00_UTC */
  char tstartstr[100], tscrstr[100];

  double tstart, tepoch, tstep, tround, tstop, tstartin, tstopin, tstepin, tstartuse, tstopuse, nseconds, chunksecs;
  char *ttotal, *tchunk;
  int ndt;
  int lmin, lmax, lminin, lmaxin, lminuse, lmaxuse, lchunk, lchunksize, lchunkfirst, lchunklast;
  int ntimechunks, nmodes, npts, imode, itime;
  int out_time_offset, out_modes_offset, out_offset, in_time_offset, in_modes_offset, in_offset, out_index, in_index;
  int iset, lminout, lmaxout;
  double tstartout, tstopout;
  float *arrptr;

  int errbufstat=setvbuf(stderr, NULL, _IONBF, BUFSIZ);
  int outbufstat=setvbuf(stdout, NULL, _IONBF, BUFSIZ);

  double wt0, wt1, wt2, wt3, wt;
  double ut0, ut1, ut2, ut3, ut;
  double st0, st1, st2, st3, st;
  double ct0, ct1, ct2, ct3, ct;

  wt0=getwalltime();
  ct0=getcputime(&ut0, &st0);

  inrecquery = (char *)cmdparams_save_str(&cmdparams, "in", &newstat);
  outseries = (char *)cmdparams_save_str(&cmdparams, "out", &newstat);
  segnamein = (char *)cmdparams_save_str(&cmdparams, "segin", &newstat);
  segnameout = (char *)cmdparams_save_str(&cmdparams, "segout", &newstat);
  int seginflag = strcmp(kNOTSPECIFIED, segnamein);
  int segoutflag = strcmp(kNOTSPECIFIED, segnameout);
  int verbflag = cmdparams_save_int(&cmdparams, "VERB", &newstat);
  int permflag = cmdparams_save_int(&cmdparams, "PERM", &newstat);
  if (permflag)
    lifetime = DRMS_PERMANENT;
  else
    lifetime = DRMS_TRANSIENT;

  char *histlinkname = (char *)cmdparams_save_str(&cmdparams, "histlink", &newstat);

  tstart=cmdparams_save_time(&cmdparams, "TSTART", &newstat);
  sprint_time(tstartstr, tstart, "TAI", 0);
  ttotal=(char *)cmdparams_save_str(&cmdparams, "TTOTAL", &newstat);
  status=drms_names_parseduration(&ttotal, &nseconds, 1);
  if (status != DRMS_SUCCESS)
  {
//    newstat = newstat | CPSAVE_UNKNOWN_ERROR;
    fprintf(stderr, "ERROR: problem parsing TTOTAL, = %s\n", ttotal);
    return 1; 
  }
  tchunk=(char *)cmdparams_save_str(&cmdparams, "TCHUNK", &newstat);
  if (strcmp(kNOTSPECIFIED, tchunk))
  {
    status=drms_names_parseduration(&tchunk, &chunksecs, 1);
    if (status != DRMS_SUCCESS)
      newstat = newstat | CPSAVE_UNKNOWN_ERROR;
  }
  else
    chunksecs=0;

  lmin=cmdparams_save_int(&cmdparams, "LMIN", &newstat);
  lmax=cmdparams_save_int(&cmdparams, "LMAX", &newstat);
  lchunksize=cmdparams_save_int(&cmdparams, "LCHUNK", &newstat);
  if (lchunksize == 0)
    lchunksize=lmax+1;

  if (newstat) 
  {
    fprintf(stderr, "ERROR: problem with input arguments, status = %d, diagnosis follows\n", newstat);
    cpsave_decode_error(newstat);
    return 1;
  }  
  else if (savestrlen != strlen(savestr)) 
  {
    fprintf(stderr, "ERROR: problem with savestr, savestrlen = %d, strlen(savestr) = %d\n", savestrlen, (int)strlen(savestr));
    return 1;
  }

  DRMS_Record_t *tempoutrec = drms_create_record(drms_env, 
                                                 outseries,
                                                 DRMS_TRANSIENT, 
                                                 &status);

  if (status != DRMS_SUCCESS) 
  {
   fprintf(stderr,"ERROR: couldn't open a record in output dataseries %s, status = %d\n", outseries, status);
   return 1;
  }

// set up ancillary dataseries for processing metadata
  char *cvsinfo = strdup("$Header: /home/cvsuser/cvsroot/JSOC/proj/globalhs/apps/jretile_manytofew.c,v 1.3 2012/09/18 16:03:27 tplarson Exp $");
  DRMS_Link_t *histlink = hcon_lookup_lower(&tempoutrec->links, histlinkname);
  if (histlink != NULL) 
  {
    histrecnum=set_history(histlink, cvsinfo);
    if (histrecnum < 0)
    {
      drms_close_record(tempoutrec, DRMS_FREE_RECORD);
      return 1;
    }
  }
  else
  {
    fprintf(stderr,"WARNING: could not find history link in output dataseries\n");
  }

// these must be present in the output dataseries and variable, not links or constants
// now done in DoIt() that calls this function
/*
  char *outchecklist[] = {"T_START", "QUALITY", "LMIN", "LMAX", "NDT"};
  DRMS_Keyword_t *outkeytest;
  int itest;
  for (itest=0; itest < ARRLENGTH(outchecklist); itest++)
  {
    outkeytest = hcon_lookup_lower(&tempoutrec->keywords, outchecklist[itest]);
    if (outkeytest == NULL || outkeytest->info->islink || outkeytest->info->recscope == 1)
    {
      fprintf(stderr, "ERROR: output keyword %s is either missing, constant, or a link\n", outchecklist[itest]);
      drms_close_record(tempoutrec, DRMS_FREE_RECORD);
      return 1;
    }
  }
*/

  tepoch=drms_getkey_time(tempoutrec, "T_START_epoch", &status);
  tstep=drms_getkey_float(tempoutrec, "T_START_step", &status);
  tround=drms_getkey_float(tempoutrec, "T_START_round", &status);
  cadence=drms_getkey_float(tempoutrec, "T_STEP", &status);
  if (fmod(tstart-tepoch,tstep) > tround/2)
  {
    sprint_time(tscrstr, tepoch, "TAI", 0);
    fprintf(stderr, "ERROR: output dataseries seems incompatible with input parameters (tstep must divide tstart-tepoch): TSTART = %s, T_START_epoch = %s, tstep = %f\n", 
                                                                                                                          tstartstr, tscrstr, tstep);
    drms_close_record(tempoutrec, DRMS_FREE_RECORD);
    return 1;
  }
  if (chunksecs == 0.0)
    chunksecs = tstep;
  else if (fmod(chunksecs,tstep))
  {
    fprintf(stderr, "ERROR: output dataseries seems incompatible with input parameters (tstep must divide chunksecs): chunksecs = %f, tstep = %f\n", chunksecs, tstep);
    drms_close_record(tempoutrec, DRMS_FREE_RECORD);
    return 1;
  }
  if (fmod(nseconds,chunksecs) != 0.0)
  {
    fprintf(stderr, "ERROR: input parameters seem incompatible (chunksecs must divide totalsecs): totalsecs = %f, chunksecs = %f\n", nseconds, chunksecs);
    drms_close_record(tempoutrec, DRMS_FREE_RECORD);
    return 1;
  }
  ntimechunks=nseconds/chunksecs;
  ndt=chunksecs/cadence;
  if (verbflag)
  {
    printf("%d timechunks, %.1f seconds per chunk\n", ntimechunks, chunksecs);
  }

  int mapmmaxout=-1;
  int sinbdivsout=-1;
  DRMS_Keyword_t *outkeytest = hcon_lookup_lower(&tempoutrec->keywords, "MAPMMAX");
  if (outkeytest != NULL && outkeytest->info->recscope == 1)
    mapmmaxout=drms_getkey_int(tempoutrec, "MAPMMAX", &status);
  outkeytest = hcon_lookup_lower(&tempoutrec->keywords, "SINBDIVS");
  if (outkeytest != NULL && outkeytest->info->recscope == 1)
    sinbdivsout=drms_getkey_int(tempoutrec, "SINBDIVS", &status);

  drms_close_record(tempoutrec, DRMS_FREE_RECORD);

  char *inchecklist[] = {"T_START", "QUALITY", "LMIN", "LMAX", "T_STEP"};
  DRMS_Keyword_t *inkeytest;
  int itest;
  inrecset = drms_open_recordset(drms_env, inrecquery, &status);
//  inrecset = drms_open_records(drms_env, inrecquery, &status);

  if (status != DRMS_SUCCESS || inrecset == NULL)
  {
    fprintf(stderr, "ERROR: problem opening input recordset: status = %d\n", status);
    return 1;
  }
  nrecsin = inrecset->n;

  if (verbflag) 
    printf("input recordset opened, nrecs = %d\n", nrecsin);

  int noinput=0;
  if (nrecsin == 0)
  {
    printf("WARNING: input recordset contains no records\n");
    noinput=1;
    goto skip1;
//    return 1;
  }

  inrec = drms_recordset_fetchnext(drms_env, inrecset, &fetchstat, &chunkstat, NULL);
//  inrec = inrecset->records[0];

  for (itest=0; itest < ARRLENGTH(inchecklist); itest++)
  {
    inkeytest = hcon_lookup_lower(&inrec->keywords, inchecklist[itest]);
    if (inkeytest == NULL)
    {
      fprintf(stderr, "ERROR: required input keyword %s is missing\n", inchecklist[itest]);
      drms_close_records(inrecset, DRMS_FREE_RECORD);
      return 1;
    }
  }

  if (cadence != drms_getkey_float(inrec, "T_STEP", &status))
  {
    fprintf(stderr, "ERROR: input T_STEP does not equal output T_STEP\n");
    drms_close_records(inrecset, DRMS_FREE_RECORD);
    return 1;
  }

  inkeytest = hcon_lookup_lower(&inrec->keywords, "MAPMMAX");
  if (inkeytest != NULL)
    mapmmax=drms_getkey_int(inrec, "MAPMMAX", &status);
  if (mapmmaxout != -1 && mapmmaxout != mapmmax)
  {
    fprintf(stderr, "ERROR: input MAPMMAX does not equal output MAPMMAX, in=%d, out=%d\n", mapmmax, mapmmaxout);
    drms_close_records(inrecset, DRMS_FREE_RECORD);
    return 1;
  }

  inkeytest = hcon_lookup_lower(&inrec->keywords, "SINBDIVS");
  if (outkeytest != NULL)
    sinbdivs=drms_getkey_int(inrec, "SINBDIVS", &status);
  if (sinbdivsout != -1 && sinbdivsout != sinbdivs)
  {
    fprintf(stderr, "ERROR: input SINBDIVS does not equal output SINBDIVS, in=%d, out=%d\n", sinbdivs, sinbdivsout);
    drms_close_records(inrecset, DRMS_FREE_RECORD);
    return 1;
  }

  status=drms_stage_records(inrecset, 1, 0);
  if (status != DRMS_SUCCESS)
  {
    fprintf(stderr, "ERROR: drms_stage_records returned status = %d\n", status);
    return 1;
  }

  skip1:

  lchunkfirst = lmin/lchunksize;
  lchunklast = lmax/lchunksize;

  nlchunks = (lchunklast - lchunkfirst) + 1;
  nrecsout = nlchunks*ntimechunks;
  outrecset = drms_create_records(drms_env, nrecsout, outseries, lifetime, &status);
  if (status != DRMS_SUCCESS || outrecset == NULL)
  {
    fprintf(stderr,"ERROR: unable to create records record in output dataseries %s, status = %d\n", outseries, status);
    drms_close_records(inrecset, DRMS_FREE_RECORD);
    return 1;
  }

/*
  irecout=0;
  for (iset=0; iset < ntimechunks; iset++)
  {
    tstartout=tstart + iset * chunksecs;
    tstopout=tstartout+chunksecs;

    for (lchunk = lchunkfirst; lchunk <= lchunklast; lchunk++)
    {
      lminout = lchunk * lchunksize; 
      lmaxout = lminout + lchunksize - 1;
      lminout = MAXIMUM(lminout,lmin);
      lmaxout = MINIMUM(lmaxout,lmax);

      outrec = outrecset->records[irecout];
      if (histlink)
        drms_setlink_static(outrec, histlinkname,  histrecnum);
      drms_setkey_int(outrec, "LMIN", lminout);
      drms_setkey_int(outrec, "LMAX", lmaxout);
      drms_setkey_time(outrec, "T_START", tstartout);
      drms_setkey_time(outrec, "T_STOP", tstopout);
      drms_setkey_time(outrec, "T_OBS", tstartout+chunksecs/2);
      drms_setkey_int(outrec, "NDT", ndt);

      irecout++;
    }
  }
*/
  long long calversout, calvers;
  int calversunset=1;

  int firsttimethrough=1;
  int *nskiparr=(int *)calloc(nrecsout,sizeof(int));
  for (irecin=0; irecin < nrecsin; irecin++)
  {
// move to end of loop when using drms_recordset_fetchnext(drms_env, inrecset, &fetchstat, &chunkstat, NULL);
//    inrec = inrecset->records[irecin]; 
    tstartin=drms_getkey_time(inrec, "T_START", &status);
    tstopin=drms_getkey_time(inrec, "T_STOP", &status);
    lminin=drms_getkey_int(inrec, "LMIN", &status);
    lmaxin=drms_getkey_int(inrec, "LMAX", &status);
    tstepin=tstopin-tstartin;

    quality=drms_getkey_int(inrec, "QUALITY", &status);
    if (status != DRMS_SUCCESS || (quality & QUAL_NODATA)) //may want stricter test on quality here
    {
      if (verbflag > 2)
      {
        sprint_time(tscrstr, tstartin, "TAI", 0);
        fprintf(stderr, "WARNING: input data not used due to quality: T_START = %s, LMIN = %d, LMAX = %d, recnum = %lld, irec = %d, status = %d, quality = %08x\n", 
                                                                      tscrstr, lminin, lmaxin, inrec->recnum, irecin, status, quality);
      }
      for (irecout=0; irecout < nrecsout; irecout++)
        nskiparr[irecout]++;
      goto continue_outer_loop;
//      continue;
    }

    if (calversunset)
    {
      calversout=drms_getkey_longlong(inrec, "CALVER64", &status);
      if (status != DRMS_SUCCESS)
        calversout = 0;
      calversunset=0;
    }

    calvers=drms_getkey_longlong(inrec, "CALVER64", &status);
    if (status != DRMS_SUCCESS)
      calvers = 0;

    if (calvers != calversout)
    {
      fprintf(stderr, "ERROR: input data has mixed CALVER64: %lld and %lld, recnum = %lld, histrecnum = %lld\n", calversout, calvers, inrec->recnum, histrecnum);
      return 0;
    }

    if (seginflag)
      segin = drms_segment_lookup(inrec, segnamein);
    else
      segin = drms_segment_lookupnum(inrec, 0);
    if (segin != NULL)
      inarr = drms_segment_read(segin, usetype, &status);
//      inarr = drms_segment_readslice(segin, usetype, startind, endind, &status);
    if (status != DRMS_SUCCESS || inarr == NULL || segin == NULL)
    {
      sprint_time(tscrstr, tstartin, "TAI", 0);
      fprintf(stderr, "ERROR: problem reading input segment, T_START = %s, LMIN = %d, LMAX = %d, recnum = %lld, irec = %d, status = %d\n", 
                                                             tscrstr, lminin, lmaxin, inrec->recnum, irecin, status);
      drms_close_records(inrecset, DRMS_FREE_RECORD);
      drms_close_records(outrecset, DRMS_FREE_RECORD);
      return 0;
    }
    else
    {
      inptr=(float *)(inarr->data);
    }

    irecout=0;
    for (iset=0; iset < ntimechunks; iset++)
    {
      tstartout=tstart + iset * chunksecs;
      tstopout=tstartout+chunksecs;
      sprint_time(tstartstr, tstartout, "TAI", 0);

      for (lchunk = lchunkfirst; lchunk <= lchunklast; lchunk++)
      {
        lminout = lchunk * lchunksize; 
        lmaxout = lminout + lchunksize - 1;
        lminout = MAXIMUM(lminout,lmin);
        lmaxout = MINIMUM(lmaxout,lmax);
        outrec = outrecset->records[irecout];

        if (firsttimethrough)
        {
          if (histlink != NULL)
            drms_setlink_static(outrec, histlinkname,  histrecnum);
          drms_copykeys(outrec, inrec, 0, kDRMS_KeyClass_Explicit);
// copykeys takes care of MAPMMAX, SINBDIVS, CALVER64, etc. that should be constant across the input
          drms_setkey_int(outrec, "LMIN", lminout);
          drms_setkey_int(outrec, "LMAX", lmaxout);
          drms_setkey_time(outrec, "T_START", tstartout);
          drms_setkey_time(outrec, "T_STOP", tstopout);
          drms_setkey_time(outrec, "T_OBS", tstartout+chunksecs/2);
          drms_setkey_int(outrec, "NDT", ndt);
        }

        if (tstartin >= tstopout || tstopin <= tstartout || lminin > lmaxout || lmaxin < lminout)
        {
          nskiparr[irecout++]++;
          continue;
        }

        if (segoutflag)
          segout = drms_segment_lookup(outrec, segnameout);
        else
          segout = drms_segment_lookupnum(outrec, 0);
        tstartuse=MAXIMUM(tstartout, tstartin);
        tstopuse= MINIMUM(tstopout, tstopin);
        lminuse=MAXIMUM(lminout, lminin);
        lmaxuse=MINIMUM(lmaxout, lmaxin);
        nmodes=MODES(lmaxuse+1)-MODES(lminuse);
        npts=(tstopuse - tstartuse)/cadence;

        out_time_offset = (tstartuse - tstartout)/cadence;
        out_modes_offset = MODES(lminuse) - MODES(lminout);
//        out_offset = 2 * (out_modes_offset * ndt + out_time_offset);
        out_offset = 0; // 2 * (out_modes_offset * npts + out_time_offset);
        in_time_offset = (tstartuse - tstartin)/cadence;
        in_modes_offset = MODES(lminuse) - MODES(lminin);
        in_offset =  2 * (in_modes_offset * tstepin / cadence + in_time_offset);

        startind[0]=2*out_time_offset;
        startind[1]=out_modes_offset;
        endind[0]=2*(out_time_offset + npts) - 1;
        endind[1]=out_modes_offset + nmodes - 1;
        totallength[0]=2*ndt;
        totallength[1]=lmaxout*(lmaxout+1)/2+lmaxout - lminout*(lminout+1)/2 + 1;

        length[0]=2*npts;
        length[1]=nmodes;
        arrptr=(float *)(calloc(length[0]*length[1],sizeof(float)));
        outarr = drms_array_create(usetype, 2, length, arrptr, &status);
        if (status != DRMS_SUCCESS || outarr == NULL || arrptr == NULL)
        {
          fprintf(stderr,"ERROR: problem creating output array: T_START = %s, LMIN = %d, LMAX = %d, length = [%d, %d], status = %d, histrecnum = %lld\n", 
                                                                tstartstr, lminout, lmaxout, length[0], length[1], status, histrecnum);
          drms_close_records(inrecset, DRMS_FREE_RECORD);
          drms_close_records(outrecset, DRMS_FREE_RECORD);
          return 0; 
        }
        outptr = (float *)(outarr->data);

        for (imode=0; imode<nmodes; imode++)
        {
          for (itime=0; itime<npts; itime++)
          {
            in_index=in_offset + 2*itime;
            out_index=out_offset + 2*itime;
            outptr[out_index] = inptr[in_index];
            outptr[out_index+1] = inptr[in_index+1];
          }
          out_offset+=2*npts;    // 2*ndt;
          in_offset+=2*tstepin/cadence;
        }

        outarr->bzero=segout->bzero;
        outarr->bscale=segout->bscale;
        status=drms_segment_writeslice_ext(segout, outarr, startind, endind, totallength, 0);
        if (status != DRMS_SUCCESS)
        {
          fprintf(stderr, "ERROR: problem writing output segment: status = %d, T_START = %s, LMIN = %d, LMAX = %d, histrecnum = %lld\n", 
                                                                  status, tstartstr, lminout, lmaxout, histrecnum);
          drms_close_records(inrecset, DRMS_FREE_RECORD);
          drms_close_records(outrecset, DRMS_FREE_RECORD);
          return 0;
        }

        drms_free_array(outarr);

        irecout++;
      }  // end loop on lchunk
    }  // end loop on iset

  firsttimethrough=0;
  drms_free_array(inarr);
  continue_outer_loop:
  inrec = drms_recordset_fetchnext(drms_env, inrecset, &fetchstat, &chunkstat, NULL);
  }  // end loop on irecin

  drms_close_records(inrecset, DRMS_FREE_RECORD);

  int nsegments=0;
  for (irecout=0; irecout < nrecsout; irecout++)
  {
    outrec=outrecset->records[irecout];
    if (noinput || nskiparr[irecout] == nrecsin)
    {
      drms_setkey_int(outrec, "QUALITY", QUAL_NODATA);
    }
    else
    {
      drms_setkey_int(outrec, "QUALITY", 0);
      nsegments++;
    }

    tnow = (double)time(NULL);
    tnow += UNIX_epoch;
    drms_setkey_time(outrec, "DATE", tnow);
  }

  free(nskiparr);
  drms_close_records(outrecset, DRMS_INSERT_RECORD);

  wt=getwalltime();
  ct=getcputime(&ut, &st);
  if (verbflag) 
  {
    printf("number of records created  = %d\n", nrecsout);
    printf("number of segments created = %d\n", nsegments);
    fprintf(stdout, "total time spent: %.2f ms wall time, %.2f ms cpu time\n", 
            wt-wt0, ct-ct0);
  }

  printf("module %s successful completion\n", cmdparams.argv[0]);

  return 0;

}

Karen Tian
Powered by
ViewCVS 0.9.4