util/apps/count_unslotted.c

/*

 count_unslotted - in=recset, counts records matching query.  Query may include the '@cadence' construct
 even in case of aia.lev1 and hmi.lev1 and other mostly empty slotted or unslotted series.

 Should be modified to simply count the lines in the input_recset file, or popen show_info -c in the general case.

*/
#include "jsoc.h"
#include "jsoc_main.h"
#include "atoinc.h"

char *module_name = "count_unslotted";

#define DIE(msg) {fprintf(stderr,"%s  Status=%d\n",msg, status); return(status?status:1);}
#define DIE_get_recset(msg) {fprintf(stderr,"$$$$ %s: %s\n", module_name, msg); return NULL;}

ModuleArgs_t module_args[] =
  {
    {ARG_STRING, "in", "NOTSPECIFIED", "input series. e.g 'mdi.fd_M_96m_lev18'"},
    {ARG_TIME, "epoch", "1993.01.01_TAI", "Reference epoch for '@' slots" },
    {ARG_END}
  };

char *get_input_recset(DRMS_Env_t *drms_env, char *inQuery);

int DoIt(void)
{
  CmdParams_t *params = &cmdparams;
  DRMS_RecordSet_t *inRS = NULL;
  int status = DRMS_SUCCESS, nrecs; 
  char inQuery[DRMS_MAXQUERYLEN];
  char in[DRMS_MAXQUERYLEN];
  const char *ingiven = params_get_str(params, "in");
  char *new_in, *in_filename;
  strcpy(in, ingiven);
  new_in = get_input_recset(drms_env, in);
  if (new_in != in)
    {
    strcpy(in, new_in);
    in_filename = new_in+1;
    }

  inRS = drms_open_records(drms_env, in, &status);
  if (status)
     {
     fprintf(stderr,"Query is: %s\n",in);
     DIE("No input data found");
     }
  nrecs = inRS->n;
  drms_close_records(inRS, DRMS_FREE_RECORD); 
  if (in_filename)
    {
    unlink(in_filename);
    }

  printf("%d\n", nrecs);
  return DRMS_SUCCESS;
}


// In cases known to not have compact slotted series and cadence is specified
// generate explicit recordset list of closest good record to desired grid
// First get vector of times and quality
// Then if vector is not OK, quit.
// then: make temp file to hold recordset list
//       start with first time to define desired grid,
//       make array of desired times.
//       make empty array of recnums
//       search vector for good images nearest desired times
//       for each found time, write record query


char *get_input_recset(DRMS_Env_t *drms_env, char *inQuery)
  {
  static char newInQuery[DRMS_MAXSERIESNAMELEN+2];
  int epoch_given = cmdparams_exists(&cmdparams, "epoch");
  TIME epoch, t_epoch;
  DRMS_Array_t *data;
  DRMS_Record_t *inTemplate;
  TIME t_start, t_stop, t_now, t_want, t_diff, this_t_diff;
  int status = 1;
  int nrecs, irec;
  int nslots, islot;
  long long *recnums;
  TIME *t_this, half;
  TIME cadence;
  double *drecnum, *dquality;
  int quality;
  long long recnum;
  char keylist[DRMS_MAXQUERYLEN];
  char filename[DRMS_MAXSERIESNAMELEN];
  char *tmpdir;
  FILE *tmpfile;
  char newIn[DRMS_MAXQUERYLEN];
  char seriesname[DRMS_MAXQUERYLEN];
  char *lbracket;
  char *at = index(inQuery, '@');
  int npkeys;
  char *timekeyname;
  double t_step;

  strcpy(seriesname, inQuery);
  lbracket = index(seriesname,'[');
  if (lbracket) *lbracket = '\0';
  inTemplate = drms_template_record(drms_env, seriesname, &status);
  if (status || !inTemplate) DIE_get_recset("Input series can not be found");

  // Now find the prime time keyword name
  npkeys = inTemplate->seriesinfo->pidx_num;
  timekeyname = NULL;
  if (npkeys > 0)
    {
    int i;
    for (i=0; i<npkeys; i++)
        {
        DRMS_Keyword_t *pkey, *skey;
        pkey = inTemplate->seriesinfo->pidx_keywords[i];
        if (pkey->info->recscope > 1)
           pkey = drms_keyword_slotfromindex(pkey);
        if (pkey->info->type != DRMS_TYPE_TIME)
          continue;
        if(i > 0) DIE_get_recset("Input series must have TIME keyword first, for now...");
        timekeyname = pkey->info->name;
        t_step = drms_keyword_getdouble(drms_keyword_stepfromslot(pkey), &status);
        if (status) DIE_get_recset("problem getting t_step");
        t_epoch = drms_keyword_getdouble(drms_keyword_epochfromslot(pkey), &status);
        if (status) DIE_get_recset("problem getting t_epoch");
        }
    }
  else
    DIE_get_recset("Must have time prime key");
  epoch = epoch_given ? params_get_time(&cmdparams, "epoch") : t_epoch;

  if (at && *at && ((strncmp(inQuery,"aia.lev1[", 9)==0 ||
                    strncmp(inQuery,"hmi.lev1[", 9)==0 ||
                    strncmp(inQuery,"aia.lev1_nrt2[",14)==0 ||
                    strncmp(inQuery,"hmi.lev1_nrt[", 13)==0 ) ||
                   epoch_given))
    {
    char *ip=(char *)inQuery, *op=newIn, *p;
    long n, mul;
    while ( *ip && ip<at )
      *op++ = *ip++;
    ip++; // skip the '@'
    n = strtol(ip, &p, 10); // get digits only
    if (*p == 's') mul = 1;
    else if (*p == 'm') mul = 60;
    else if (*p == 'h') mul = 3600;
    else if (*p == 'd') mul = 86400;
    else 
      DIE_get_recset("cant make sense of @xx cadence spec");
    cadence = n * mul;
    ip = ++p;  // skip cadence multiplier
    while ( *ip )
      *op++ = *ip++;
    *op = '\0';
    half = cadence/2.0;
    sprintf(keylist, "%s,QUALITY,recnum", timekeyname);
    data = drms_record_getvector(drms_env, newIn, keylist, DRMS_TYPE_DOUBLE, 0, &status);
    if (!data || status)
      {
      fprintf(stderr, "status=%d\n", status);
      DIE_get_recset("getkey_vector failed\n");
      }
    nrecs = data->axis[1];
    irec = 0;
    t_this = (TIME *)data->data;
    dquality = (double *)data->data + 1*nrecs;
    drecnum = (double *)data->data + 2*nrecs;
    if (epoch_given)
      {
      int s0 = (t_this[0] - epoch)/cadence;
      TIME t0 = s0*cadence + epoch;
      t_start = (t0 < t_this[0] ? t0 + cadence : t0);
      }
    else
      t_start = t_this[0];
    t_stop = t_this[nrecs-1];
    nslots = (t_stop - t_start + cadence/2)/cadence;
    recnums = (long long *)malloc(nslots*sizeof(long long));
    for (islot=0; islot<nslots; islot++)
      recnums[islot] = 0;
    islot = 0;
    t_want = t_start;
    t_diff = 1.0e9;
    for (irec = 0; irec<nrecs; irec++)
        {
        t_now = t_this[irec];
        quality = (int)dquality[irec] & 0xFFFFFFFF;
        recnum = (long long)drecnum[irec];
        this_t_diff = fabs(t_now - t_want);
        if (quality < 0)
          continue;
        if (t_now <= (t_want-half))
          continue;
        while (t_now > (t_want+half))
          {
          islot++;
          if (islot >= nslots)
             break;
          t_want = t_start + cadence * islot;
          this_t_diff = fabs(t_now - t_want);
          t_diff = 1.0e8;
          }
        if (islot < nslots && this_t_diff <= t_diff)
          recnums[islot] = recnum;
        t_diff = fabs(t_now - t_want);
        }
    if (islot+1 < nslots)
      nslots = islot+1;  // take what we got.
    tmpdir = getenv("TMPDIR");
    if (!tmpdir) tmpdir = "/tmp";
    sprintf(filename, "%s/%sXXXXXX", tmpdir, module_name);
    mkstemp(filename);
    tmpfile = fopen(filename,"w");
    for (islot=0; islot<nslots; islot++)
      if (recnums[islot])
        fprintf(tmpfile, "%s[:#%lld]\n", seriesname, recnums[islot]);
    fclose(tmpfile);
    free(recnums);
    drms_free_array(data);
    sprintf(newInQuery,"@%s", filename);
    return(newInQuery);
    }
  else
          return(inQuery);
  }
Revision:	1.1
Committed:	Mon Apr 25 01:09:36 2016 UTC (7 years, 5 months ago) by phil
Content type:	text/plain
Branch:	MAIN
CVS Tags:	Ver_8-11, Ver_LATEST, Ver_9-41, Ver_8-12, Ver_9-5, Ver_9-4, Ver_9-3, Ver_9-2, Ver_9-1, Ver_9-0, HEAD
Log Message:	Initial add. Counts records in the get_input_recset call for cases with mostly unfilled slots. Useful for aia.lev1 and hmi.lev1 where a reduced cadence is being queried.
#	Content
1	/*
2
3	count_unslotted - in=recset, counts records matching query. Query may include the '@cadence' construct
4	even in case of aia.lev1 and hmi.lev1 and other mostly empty slotted or unslotted series.
5
6	Should be modified to simply count the lines in the input_recset file, or popen show_info -c in the general case.
7
8	*/
9	#include "jsoc.h"
10	#include "jsoc_main.h"
11	#include "atoinc.h"
12
13	char *module_name = "count_unslotted";
14
15	#define DIE(msg) {fprintf(stderr,"%s Status=%d\n",msg, status); return(status?status:1);}
16	#define DIE_get_recset(msg) {fprintf(stderr,"$$$$ %s: %s\n", module_name, msg); return NULL;}
17
18	ModuleArgs_t module_args[] =
19	{
20	{ARG_STRING, "in", "NOTSPECIFIED", "input series. e.g 'mdi.fd_M_96m_lev18'"},
21	{ARG_TIME, "epoch", "1993.01.01_TAI", "Reference epoch for '@' slots" },
22	{ARG_END}
23	};
24
25	char get_input_recset(DRMS_Env_t drms_env, char *inQuery);
26
27	int DoIt(void)
28	{
29	CmdParams_t *params = &cmdparams;
30	DRMS_RecordSet_t *inRS = NULL;
31	int status = DRMS_SUCCESS, nrecs;
32	char inQuery[DRMS_MAXQUERYLEN];
33	char in[DRMS_MAXQUERYLEN];
34	const char *ingiven = params_get_str(params, "in");
35	char new_in, in_filename;
36	strcpy(in, ingiven);
37	new_in = get_input_recset(drms_env, in);
38	if (new_in != in)
39	{
40	strcpy(in, new_in);
41	in_filename = new_in+1;
42	}
43
44	inRS = drms_open_records(drms_env, in, &status);
45	if (status)
46	{
47	fprintf(stderr,"Query is: %s\n",in);
48	DIE("No input data found");
49	}
50	nrecs = inRS->n;
51	drms_close_records(inRS, DRMS_FREE_RECORD);
52	if (in_filename)
53	{
54	unlink(in_filename);
55	}
56
57	printf("%d\n", nrecs);
58	return DRMS_SUCCESS;
59	}
60
61
62	// In cases known to not have compact slotted series and cadence is specified
63	// generate explicit recordset list of closest good record to desired grid
64	// First get vector of times and quality
65	// Then if vector is not OK, quit.
66	// then: make temp file to hold recordset list
67	// start with first time to define desired grid,
68	// make array of desired times.
69	// make empty array of recnums
70	// search vector for good images nearest desired times
71	// for each found time, write record query
72
73
74	char get_input_recset(DRMS_Env_t drms_env, char *inQuery)
75	{
76	static char newInQuery[DRMS_MAXSERIESNAMELEN+2];
77	int epoch_given = cmdparams_exists(&cmdparams, "epoch");
78	TIME epoch, t_epoch;
79	DRMS_Array_t *data;
80	DRMS_Record_t *inTemplate;
81	TIME t_start, t_stop, t_now, t_want, t_diff, this_t_diff;
82	int status = 1;
83	int nrecs, irec;
84	int nslots, islot;
85	long long *recnums;
86	TIME *t_this, half;
87	TIME cadence;
88	double drecnum, dquality;
89	int quality;
90	long long recnum;
91	char keylist[DRMS_MAXQUERYLEN];
92	char filename[DRMS_MAXSERIESNAMELEN];
93	char *tmpdir;
94	FILE *tmpfile;
95	char newIn[DRMS_MAXQUERYLEN];
96	char seriesname[DRMS_MAXQUERYLEN];
97	char *lbracket;
98	char *at = index(inQuery, '@');
99	int npkeys;
100	char *timekeyname;
101	double t_step;
102
103	strcpy(seriesname, inQuery);
104	lbracket = index(seriesname,'[');
105	if (lbracket) *lbracket = '\0';
106	inTemplate = drms_template_record(drms_env, seriesname, &status);
107	if (status \|\| !inTemplate) DIE_get_recset("Input series can not be found");
108
109	// Now find the prime time keyword name
110	npkeys = inTemplate->seriesinfo->pidx_num;
111	timekeyname = NULL;
112	if (npkeys > 0)
113	{
114	int i;
115	for (i=0; i<npkeys; i++)
116	{
117	DRMS_Keyword_t pkey, skey;
118	pkey = inTemplate->seriesinfo->pidx_keywords[i];
119	if (pkey->info->recscope > 1)
120	pkey = drms_keyword_slotfromindex(pkey);
121	if (pkey->info->type != DRMS_TYPE_TIME)
122	continue;
123	if(i > 0) DIE_get_recset("Input series must have TIME keyword first, for now...");
124	timekeyname = pkey->info->name;
125	t_step = drms_keyword_getdouble(drms_keyword_stepfromslot(pkey), &status);
126	if (status) DIE_get_recset("problem getting t_step");
127	t_epoch = drms_keyword_getdouble(drms_keyword_epochfromslot(pkey), &status);
128	if (status) DIE_get_recset("problem getting t_epoch");
129	}
130	}
131	else
132	DIE_get_recset("Must have time prime key");
133	epoch = epoch_given ? params_get_time(&cmdparams, "epoch") : t_epoch;
134
135	if (at && *at && ((strncmp(inQuery,"aia.lev1[", 9)==0 \|\|
136	strncmp(inQuery,"hmi.lev1[", 9)==0 \|\|
137	strncmp(inQuery,"aia.lev1_nrt2[",14)==0 \|\|
138	strncmp(inQuery,"hmi.lev1_nrt[", 13)==0 ) \|\|
139	epoch_given))
140	{
141	char ip=(char )inQuery, op=newIn, p;
142	long n, mul;
143	while ( *ip && ip<at )
144	op++ = ip++;
145	ip++; // skip the '@'
146	n = strtol(ip, &p, 10); // get digits only
147	if (*p == 's') mul = 1;
148	else if (*p == 'm') mul = 60;
149	else if (*p == 'h') mul = 3600;
150	else if (*p == 'd') mul = 86400;
151	else
152	DIE_get_recset("cant make sense of @xx cadence spec");
153	cadence = n * mul;
154	ip = ++p; // skip cadence multiplier
155	while ( *ip )
156	op++ = ip++;
157	*op = '\0';
158	half = cadence/2.0;
159	sprintf(keylist, "%s,QUALITY,recnum", timekeyname);
160	data = drms_record_getvector(drms_env, newIn, keylist, DRMS_TYPE_DOUBLE, 0, &status);
161	if (!data \|\| status)
162	{
163	fprintf(stderr, "status=%d\n", status);
164	DIE_get_recset("getkey_vector failed\n");
165	}
166	nrecs = data->axis[1];
167	irec = 0;
168	t_this = (TIME *)data->data;
169	dquality = (double )data->data + 1nrecs;
170	drecnum = (double )data->data + 2nrecs;
171	if (epoch_given)
172	{
173	int s0 = (t_this[0] - epoch)/cadence;
174	TIME t0 = s0*cadence + epoch;
175	t_start = (t0 < t_this[0] ? t0 + cadence : t0);
176	}
177	else
178	t_start = t_this[0];
179	t_stop = t_this[nrecs-1];
180	nslots = (t_stop - t_start + cadence/2)/cadence;
181	recnums = (long long )malloc(nslotssizeof(long long));
182	for (islot=0; islot<nslots; islot++)
183	recnums[islot] = 0;
184	islot = 0;
185	t_want = t_start;
186	t_diff = 1.0e9;
187	for (irec = 0; irec<nrecs; irec++)
188	{
189	t_now = t_this[irec];
190	quality = (int)dquality[irec] & 0xFFFFFFFF;
191	recnum = (long long)drecnum[irec];
192	this_t_diff = fabs(t_now - t_want);
193	if (quality < 0)
194	continue;
195	if (t_now <= (t_want-half))
196	continue;
197	while (t_now > (t_want+half))
198	{
199	islot++;
200	if (islot >= nslots)
201	break;
202	t_want = t_start + cadence * islot;
203	this_t_diff = fabs(t_now - t_want);
204	t_diff = 1.0e8;
205	}
206	if (islot < nslots && this_t_diff <= t_diff)
207	recnums[islot] = recnum;
208	t_diff = fabs(t_now - t_want);
209	}
210	if (islot+1 < nslots)
211	nslots = islot+1; // take what we got.
212	tmpdir = getenv("TMPDIR");
213	if (!tmpdir) tmpdir = "/tmp";
214	sprintf(filename, "%s/%sXXXXXX", tmpdir, module_name);
215	mkstemp(filename);
216	tmpfile = fopen(filename,"w");
217	for (islot=0; islot<nslots; islot++)
218	if (recnums[islot])
219	fprintf(tmpfile, "%s[:#%lld]\n", seriesname, recnums[islot]);
220	fclose(tmpfile);
221	free(recnums);
222	drms_free_array(data);
223	sprintf(newInQuery,"@%s", filename);
224	return(newInQuery);
225	}
226	else
227	return(inQuery);
228	}