source: CMIP6dreq/trunk/dreqPy/scope.py @ 880

Subversion URL: http://proj.badc.rl.ac.uk/svn/exarch/CMIP6dreq/trunk/dreqPy/scope.py@880
Revision 880, 59.0 KB checked in by mjuckes, 3 years ago (diff)

Updated setup for tag 01.beta.43

Line 
1"""Date Request Scoping module
2------------------------------
3The scope.py module contains the dreqQuery class and a set of ancilliary functions. The dreqQuery class contains methods for analysing the data request.
4"""
5
6try:
7  import dreq
8  imm=1
9except:
10  import dreqPy.dreq  as dreq
11  imm=2
12
13if imm == 1:
14  from utilities import cmvFilter
15  import makeTables
16  import fgrid
17  import volsum
18else:
19  import dreqPy.volsum as volsum
20  import dreqPy.fgrid as fgrid
21  from dreqPy.utilities import cmvFilter
22  import dreqPy.makeTables as makeTables
23
24import collections, string, operator
25import sys, os
26
27def sortTimeSlice( tsl ):
28  s = set()
29  for ts in tsl:
30    s.add( ts[1] )
31  if len(s) > 1:
32    return (-1,None,'Multiple slice types')
33  tst = s.pop()
34  if not ( tst == 'simpleRange' or (len(tst) > 13 and tst[:13] == 'branchedYears') ):
35    return (-2,None,'slice type aggregation not supported')
36  if len(tsl) == 2:
37    tsll = list( tsl )
38    sa,ea = tsll[0][2:]
39    sb,eb = tsll[1][2:]
40    if sa <= sb and ea >= eb:
41      return (1,tsll[0], 'Taking largest slice')
42    if sb <= sa and eb >= ea:
43      return (1,tsll[1], 'Taking largest slice')
44    if ea < sb or eb < sa:
45      return (2,tsll, 'Slices are disjoint')
46    return (-3,None, 'Overlapping slices')
47  else:
48    tsll = sorted( list(tsl), key=lambda x: x[3] )
49    if min( [x[2] for x in tsll] ) == tsll[-1][2]:
50        return (1,tsll[-1], 'Taking largest slice')
51    return (-4,None, 'Cannot sort slices')
52
53odsz = {'landUse':(5,'free'), 'tau':7, 'scatratio':15, 'effectRadLi|tau':(28,'query pending'), 'vegtype':(8,'free'), 'sza5':5, 'site':(119,'73 for aquaplanet .. '), 'iceband':(5,'free'), 'dbze':15, 'spectband':(10,'free'), 'misrBands':(7,'query pending'), 'effectRadIc|tau':(28,'query pending')}
54
55python2 = True
56if sys.version_info[0] == 3:
57  python2 = False
58  from functools import reduce
59  try: 
60    from utilP3 import mlog3
61  except:
62    from dreqPy.utilP3 import mlog3
63  mlg = mlog3()
64else:
65  from utilP2 import util
66  mlg = util.mlog()
67
68class c1(object):
69  def __init__(self):
70    self.a = collections.defaultdict( int )
71class c1s(object):
72  def __init__(self):
73    self.a = collections.defaultdict( set )
74
75NT_txtopts = collections.namedtuple( 'txtopts', ['mode'] )
76
77def vfmt(ss):
78  stb = ss*1.e-12
79  if stb < .099:
80    return '%7.2fGb' % (stb*100)
81  else:
82    return '%7.2fTb' % stb
83
84class baseException(Exception):
85  """Basic exception for general use in code."""
86
87  def __init__(self,msg):
88    self.msg = 'scope:: %s' % msg
89
90  def __str__(self):
91    return repr( self.msg )
92
93  def __repr__(self):
94    return self.msg
95
96nt_mcfg = collections.namedtuple( 'mcfg', ['nho','nlo','nha','nla','nlas','nls','nh1'] )
97class cmpd(object):
98  def __init__(self,dct):
99    self.d = dct
100  def cmp(self,x,y,):
101    return cmp( self.d[x], self.d[y] )
102
103
104def filter1( a, b ):
105  if b < 0:
106    return a
107  else:
108    return min( [a,b] )
109
110def filter2( a, b, tt, tm ):
111## largest tier less than or equal to tm
112  ll = [t for t in tt if t <= tm]
113  if len( ll ) > 0:
114    t1 = [t for t in tt if t <= tm][-1]
115    it1 = tt.index(t1)
116    aa = a[it1]
117    if b < 0:
118      return aa
119    else:
120      return min( [aa,b] )
121  else:
122    return 0
123
124npy = {'1hrClimMon':24*12, 'daily':365, u'Annual':1, u'fx':0.01, u'1hr':24*365, u'3hr':8*365,
125       u'monClim':12, u'Timestep':100, u'6hr':4*365, u'day':365, u'1day':365, u'mon':12, u'yr':1,
126       u'1mon':12, 'month':12, 'year':1, 'monthly':12, 'hr':24*365, 'other':24*365,
127        'subhr':24*365, 'Day':365, '6h':4*365, '3 hourly':8*365, '':1 }
128
129## There are 4 cmor variables with blank frequency ....
130
131def vol01( sz, v, npy, freq, inx ):
132  n1 = npy[freq]
133  s = sz[inx.uid[v].stid]
134  assert type(s) == type(1), 'Non-integer size found for %s' % v
135  assert type(n1) in (type(1),type(0.)), 'Non-number "npy" found for %s, [%s]' % (v,freq)
136  return s*n1
137
138class col_list(object):
139  def __init__(self):
140    self.a = collections.defaultdict(list)
141
142class col_count(object):
143  def __init__(self):
144    self.a = collections.defaultdict(int)
145
146class dreqQuery(object):
147  __doc__ = """Methods to analyse the data request, including data volume estimates"""
148  def __init__(self,dq=None,tierMax=1):
149    if dq == None:
150      self.dq = dreq.loadDreq()
151    else:
152      self.dq=dq
153    self.rlu = {}
154    for i in self.dq.coll['objective'].items:
155      k = '%s.%s' % (i.mip,i.label)
156      assert not k in self.rlu, 'Duplicate label in objectives: %s' % k
157      self.rlu[k] = i.uid
158
159    self.odsz = odsz
160    self.npy = npy
161    self.strSz = dict()
162    self.cmvFilter = cmvFilter( self )
163    self.tierMax = tierMax
164    self.gridPolicyDefaultNative = False
165    self.gridPolicyTopOnly = True
166    self.exptFilter = None
167    self.exptFilterBlack = None
168    self.uniqueRequest = False
169
170    self.mips = set( [x.label for x in self.dq.coll['mip'].items ] )
171    self.mips = ['CMIP','AerChemMIP', 'C4MIP', 'CFMIP', 'DAMIP', 'DCPP', 'FAFMIP', 'GeoMIP', 'GMMIP', 'HighResMIP', 'ISMIP6', 'LS3MIP', 'LUMIP', 'OMIP', 'PMIP', 'RFMIP', 'ScenarioMIP', 'VolMIP', 'CORDEX', 'DynVar', 'SIMIP', 'VIACSAB']
172    self.mipsp = self.mips[:-4]
173    self.cmvGridId, i4 = fgrid.fgrid( self.dq )
174    assert len(i4) == 0
175
176    self.experiments = set( [x.uid for x in self.dq.coll['experiment'].items ] )
177    self.exptByLabel = {}
178    for x in self.dq.coll['experiment'].items:
179      if x.label in self.exptByLabel:
180        print ( 'ERROR: experiment label duplicated: %s' % x.label )
181      self.exptByLabel[x.label] = x.uid
182    self.mipls = sorted( list( self.mips ) )
183
184    self.default_mcfg = nt_mcfg._make( [259200,60,64800,40,20,5,100] )
185    self.mcfg = self.default_mcfg._asdict()
186    self.mcfgNote = None
187    self.szcfg()
188    self.requestItemExpAll(  )
189
190  def setMcfg(self, ll, msg=None ):
191    assert len(ll) == 7, 'Model config must be of length 7: %s' % str(ll)
192    assert all( [type(x) == type(1) for x in ll] )
193    self.mcfg = nt_mcfg._make( ll )._asdict()
194    if msg == None:
195      self.mcfgNote = 'User supplied model configuration: %s' % str(ll)
196    else:
197      self.mcfgNote = msg
198    self.szcfg()
199
200  def szcfg(self):
201    szr = {'100km':64800, '1deg':64800, '2deg':16200 }
202    self.szss = {}
203    self.sz = {}
204    self.szg = collections.defaultdict( dict )
205    self.szgss = collections.defaultdict( dict )
206    self.isLatLon = {}
207    self.szSrf = collections.defaultdict( dict )
208    self.szssSrf = collections.defaultdict( dict )
209    for i in self.dq.coll['spatialShape'].items:
210      gtype = 'a'
211      if i.levelFlag == False:
212        ds =  i.dimensions.split( '|' )
213        if ds[-1] in ['site', 'basin']:
214          vd = ds[-2]
215        else:
216          vd = ds[-1]
217 
218        if vd[:4] == 'olev' or vd == 'rho':
219          gtype = 'o'
220          nz = self.mcfg['nlo']
221        elif vd[:4] == 'alev':
222          nz = self.mcfg['nla']
223        elif vd in ['slevel']:
224          nz = self.mcfg['nls']
225        elif vd in ['snowdepth','sdepth']:
226          nz = 5
227        elif vd == 'aslevel':
228          nz = self.mcfg['nlas']
229        else:
230          mlg.prnt( 'Failed to parse dimensions %s' % i.dimensions )
231          raise
232      else:
233        nz = i.levels
234
235      dims = set( i.dimensions.split( '|' ) )
236      if 'latitude' in dims and 'longitude' in dims:
237        if gtype == 'o':
238          nh = self.mcfg['nho']
239          self.isLatLon[i.uid] = 'o'
240        else:
241          nh = self.mcfg['nha']
242          self.isLatLon[i.uid] = 'a'
243      else:
244        nh = 10
245        self.isLatLon[i.uid] = False
246
247      self.szss[i.uid] = nh*nz
248      if self.isLatLon[i.uid] != False and len(dims) == 2:
249        self.szssSrf[i.uid] = { 'a':self.mcfg['nha']*nz, 'o':self.mcfg['nho']*nz }
250
251      for k in szr:
252        if self.isLatLon[i.uid] != False:
253          self.szgss[k][i.uid] = szr[k]*nz
254        else:
255          self.szgss[k][i.uid] = nh*nz
256
257    for i in self.dq.coll['structure'].items:
258      s = 1
259      knownAtmos = False
260      if i.odims != '':
261        if i.odims in odsz:
262           sf = odsz[i.odims]
263        else:
264           ## print 'SEVERE.odims.00001: no information on dimension size: %s' % i.odims
265           sf = 5
266        if type( sf ) == type( () ):
267          sf = sf[0]
268        s = s*sf
269        if i.odims not in ['iceband']:
270          knownAtmos = True
271      if i.spid in self.szss:
272        self.sz[i.uid] = self.szss[i.spid]*s
273
274        if i.uid in self.szssSrf:
275          if knownAtmos:
276            self.sz[i.uid] = self.szssSrf[i.spid]['a']*s
277          else:
278            for k in ['a','o']:
279               self.szSrf[i.uid][k] = self.szssSrf[i.spid][k]*s
280
281        for k in szr:
282          self.szg[k][i.uid] = self.szgss[k][i.spid]*s
283      else:
284        print ('WARNING: spid has no size info: %s [%s]' % (i.spid,i.uid) )
285        self.sz[i.uid] = 0.
286        for k in szr:
287          self.szg[k][i.uid] = 0.
288
289  def getRequestLinkByMip( self, mipSel ):
290    """Return the set of request links which are associated with specified MIP"""
291
292    if type(mipSel) == type( {} ):
293      return self.getRequestLinkByMipObjective(self,mipSel)
294
295    if type(mipSel) == type(''):
296      t1 = lambda x: x == mipSel
297    elif type(mipSel) == type(set()):
298      t1 = lambda x: x in mipSel
299
300    s = set()
301    for i in self.dq.coll['requestLink'].items:
302      if t1(i.mip):
303        if 'requestItem' in self.dq.inx.iref_by_sect[i.uid].a:
304          if any( [ self.rqiExp[x][3] > 0 for  x in self.dq.inx.iref_by_sect[i.uid].a['requestItem'] if x in self.rqiExp ] ):
305            s.add( i )
306
307    self.rqs = list( s )
308    return self.rqs
309
310  def getRequestLinkByMipObjective( self, mipSel ):
311    """Return the set of request links which are associated with specified MIP and its objectives"""
312
313    assert type(mipSel) == type( {} ),'Argument must be a dictionary, listing objectives for each MIP'
314
315    s = set()
316    for i in self.dq.coll['requestLink'].items:
317      if i.mip in mipSel:
318        if len(mipSel[i.mip]) == 0:
319          s.add( i )
320        elif 'objectiveLink' in self.dq.inx.iref_by_sect[i.uid].a:
321          ss = set( [self.dq.inx.uid[k].label for k in self.dq.inx.iref_by_sect[i.uid].a['objectiveLink'] ] )
322          if any( [x in mipSel[i.mip] for x in ss] ):
323            s.add( i )
324##
325## filter requestLinks by tierMax: check to see whether they link to experiments with tier below or equal to tiermax.
326##
327    s1 = set()
328    for i in s:
329      if 'requestItem' in self.dq.inx.iref_by_sect[i.uid].a:
330        if any( [ self.rqiExp[x][-1] > 0 for  x in self.dq.inx.iref_by_sect[i.uid].a['requestItem'] if x in self.rqiExp ] ):
331            s1.add( i )
332
333    self.rqs = list( s1 )
334    return self.rqs
335
336  def varGroupXexpt(self, rqList ):
337    """For a list of request links, return a list of variable group IDs for each experiment"""
338    self.cc = collections.defaultdict( list )
339    ## dummy = {self.cc[i.expt].append(i.rlid) for i in self.dq.coll['requestItem'].items if i.rlid in {j.uid for j in rqList} }
340    return self.cc
341
342  def yearsInRequest(self, rql ):
343    self.ntot = sum( [i.ny for i in self.dq.coll['requestItem'].items if i.rlid == rql.uid] )
344    return self.ntot
345
346  def rqlByExpt( self, l1, ex, pmax=2, expFullEx=False ):
347    """rqlByExpt: return a set of request links for an experiment"""
348##
349    inx = self.dq.inx
350
351    if ex != None:
352   
353      exi = self.dq.inx.uid[ex]
354      if exi._h.label == 'experiment':
355        exset = set( [ex,exi.egid,exi.mip] )
356      else:
357        exset = set( self.esid_to_exptList(ex,deref=False,full=expFullEx) )
358##
359## rql is the set of all request links which are associated with a request item for this experiment set
360##
361   
362      l1p = set()
363      for i in l1:
364        if i.preset < 0 or i.preset <= pmax:
365          if i.esid in exset:
366            l1p.add(i)
367    else:
368      exset = None
369      l1p = l1
370
371    rql0 = set()
372    for i in l1p:
373       rql0.add(i.rlid)
374
375    rqlInv = set()
376    for u in rql0:
377      if inx.uid[u]._h.label == 'remarks':
378        rqlInv.add( u )
379    if len(rqlInv) != 0:
380      mlg.prnt ( 'WARNING.001.00002: %s invalid request links from request items ...' % len(rqlInv) )
381    rql = set()
382    for u in rql0:
383       if inx.uid[u]._h.label != 'remarks':
384         rql.add( u ) 
385
386    return rql, l1p, exset
387
388  def varsByRql( self, rql, pmax=2, intersection=False, asDict=False): 
389      """The complete set of variables associated with a set of request links."""
390      inx = self.dq.inx
391      cc1 = collections.defaultdict( set )
392      for i in rql:
393        o = inx.uid[i]
394        if o.opt == 'priority':
395          p = int( float( o.opar ) )
396          assert p in [1,2,3], 'Priority incorrectly set .. %s, %s, %s' % (o.label,o.title, o.uid)
397          cc1[inx.uid[i].mip].add( (inx.uid[i].refid,p) )
398        else:
399          cc1[inx.uid[i].mip].add( inx.uid[i].refid )
400
401      if intersection:
402        ccv = {}
403#
404# set of request variables for each MIP
405#
406        for k in cc1:
407          thisc = reduce( operator.or_, [set( inx.iref_by_sect[vg].a['requestVar'] ) for vg in cc1[k] ] )
408          rqvgs = collections.defaultdict( set )
409          for x in cc1[k]:
410            if type(x) == type( () ):
411              rqvgs[x[0]].add( x[1] )
412            else:
413              rqvgs[x].add( 3 )
414         
415          s = set()
416          for vg in rqvgs:
417            for l in inx.iref_by_sect[vg].a['requestVar']:
418              if inx.uid[l].priority <= min(pmax,max(rqvgs[vg])):
419                s.add( inx.uid[l].vid )
420          ccv[k] = s
421
422        if len( ccv.keys() ) < len( list(imips) ):
423          vars = set()
424        else:
425          vars =  reduce( operator.and_, [ccv[k] for k in ccv] )
426      else:
427        rqvgs = collections.defaultdict( set )
428        for k in cc1:
429          for x in cc1[k]:
430            if type(x) == type( () ):
431              rqvgs[x[0]].add( x[1] )
432            else:
433              rqvgs[x].add( 3 )
434         
435###To obtain a set of variables associated with this collection of variable groups:
436
437        if asDict:
438          vars = collections.defaultdict( list )
439        else:
440          vars = set()
441        for vg in rqvgs:
442          for l in inx.iref_by_sect[vg].a['requestVar']:
443            if inx.uid[l].priority <= min(pmax,max(rqvgs[vg])):
444               if asDict:
445                 vars[inx.uid[l].vid].append( vg )
446               else:
447                 vars.add(inx.uid[l].vid)
448
449        ##col1 = reduce( operator.or_, [set( inx.iref_by_sect[vg].a['requestVar'] ) for vg in rqvg ] )
450### filter out cases where the request does not point to a CMOR variable.
451    ##vars = {vid for vid in vars if inx.uid[vid][0] == u'CMORvar'}
452
453      if asDict:
454        thisvars = {}
455        for vid in vars:
456           if inx.uid[vid]._h.label == u'CMORvar':
457             thisvars[vid] = vars[vid]
458      else:
459        thisvars = set()
460        for vid in vars:
461           if inx.uid[vid]._h.label == u'CMORvar':
462             thisvars.add(vid)
463
464      return thisvars
465
466  def exptYears( self, rqll, ex=None, exBlack=None):
467    """Parse a set of request links, and get years requested for each (varGroup, expt, grid) tuple """
468     
469    self.tsliceDict = collections.defaultdict( dict )
470    ccts = collections.defaultdict( dict )
471    cc = collections.defaultdict( set )
472    for rl in rqll:
473      if 'requestItem' not in self.dq.inx.iref_by_sect[rl.uid].a:
474        print ( 'WARN.001.00001: no request items for: %s, %s' % (rl.uid, rl.title) )
475      else:
476
477        if rl.grid == '100km':
478            grd = '1deg'
479        if rl.grid in ['1deg','2deg']:
480            grd = rl.grid
481        else:
482          ## note that naming of "gridreq" is unfortunate ... "No" means that native grid is required
483          if rl.gridreq in ['No', 'no'] or self.gridPolicyDefaultNative:
484            grd = 'native'
485          else:
486            ## print ( 'INFO.grd.00001: defaulting to grid ..%s, %s, %s' % (rl.label,rl.title, rl.uid) )
487            grd = 'DEF'
488
489        for iu in self.dq.inx.iref_by_sect[rl.uid].a['requestItem']:
490          i = self.dq.inx.uid[iu]
491
492##
493## apply "treset" filter to request items linked to this group.
494##
495          if self.tierMax < 0 or 'treset' not in i.__dict__ or i.treset <= self.tierMax:
496            if iu in self.rqiExp:
497              for e in self.rqiExp[iu][1]:
498                if (ex == None or e in ex) and (exBlack == None or e not in exBlack):
499                  this = self.rqiExp[iu][1][e]
500                  if this != None:
501                    thisns = this[-3]
502                    thisny = this[-2]
503                    thisne = this[-1]
504                    ##cc[ (rl.refid,e,grd) ].add( filter1( thisns*thisny*thisne, i.nymax) )
505                    cc[ (rl.refid,e,grd) ].add( thisns*thisny*thisne )
506                    if self.rqiExp[iu][4] != None:
507                      ccts[(rl.refid,e)][thisns*thisny*thisne] = self.rqiExp[iu][4]
508
509    ee = collections.defaultdict( dict )
510
511    revertToLast = True
512    if revertToLast:
513      for g,e,grd in cc:
514        ee[g][(e,grd)] = max( cc[( g,e,grd) ] )
515        if (g,e) in ccts and ee[g][(e,grd)] in ccts[(g,e)]:
516           self.tsliceDict[g][e] = ccts[(g,e)][ ee[g][(e,grd)] ]
517      return ee
518    ff = collections.defaultdict( dict )
519##
520## this needs to be done separately for ocean and atmosphere, because of the default logic workflow ...
521    for g,e,grd in cc:
522      ee[(g,e)][grd] = max( cc[( g,e,grd) ] )
523
524    xx = collections.defaultdict( dict )
525    for g,e in ee:
526      ddef = ee[(g,e)].get( 'DEF', 0 )
527      for grd in ee[(g,e)]:
528        if grd != 'DEF':
529          xx[(g,'a')][(e, grd)] = ee[(g,e)][grd]
530          xx[(g,'o')][(e, grd)] = ee[(g,e)][grd]
531          xx[(g,'')][(e, grd)] = ee[(g,e)][grd]
532        if grd == 'native' and ddef != 0:
533          xx[(g,'a')][(e, 'native')] = max( [xx[(g,'a')][(e, 'native')],ddef] )
534          xx[(g,'')][(e, 'native')] = max( [xx[(g,'')][(e, 'native')],ddef] )
535        if grd == '1deg' and ddef != 0:
536          xx[(g,'o')][(e, '1deg')] = max( [xx[(g,'o')][(e, '1deg')],ddef] )
537
538    for grp,flg in xx:
539      ff[grp][flg] = xx[(grp,flg)]
540         
541    ## return dict[<variable group>]{dict[<experiment><grid>]{<years>}}
542    ## return dict[<variable group>][grid flag]{dict[<experiment>,<grid>]{<years>}}
543    return ff
544
545  def volByExpt( self, l1, ex, pmax=1, cc=None, retainRedundantRank=False, intersection=False,expFullEx=False, adsCount=False ):
546    """volByExpt: calculates the total data volume associated with an experiment/experiment group and a list of request items.
547          The calculation has some approximations concerning the number of years in each experiment group.
548          cc: an optional collector, to accumulate indexed volumes. """
549##
550    inx = self.dq.inx
551    imips = set()
552    for i in l1:
553      imips.add(i.mip)
554   
555    rql, l1p, exset = self.rqlByExpt( l1, ex, pmax=pmax, expFullEx=expFullEx )
556    verbose = False
557    if verbose:
558      for i in rql:
559        r = inx.uid[i]
560        print ( '%s, %s, %s' % (r.label, r.title, r.uid) )
561
562    if ex != None:
563     
564      exi = self.dq.inx.uid[ex]
565      if exi._h.label == 'experiment':
566        exset = set( [ex,exi.egid,exi.mip] )
567#####
568    if len( rql ) == 0:
569      self.vars = set()
570      return (0,{},{} )
571
572## The complete set of variables associated with these requests:
573    vars = self.varsByRql( rql, pmax=pmax, intersection=intersection, asDict=True) 
574##
575## filter by configuration option and rank
576##
577    if not retainRedundantRank:
578      len1 = len(vars.keys())
579      cmv = self.cmvFilter.filterByChoiceRank(cmv=vars.keys())
580      vars = cmv
581   
582    self.vars = vars
583
584    e = {}
585    for u in rql:
586### for request variables which reference the variable group attached to the link, add the associate CMOR variables, subject to priority
587      i = inx.uid[u]
588      e[i.uid] = set()
589      si = collections.defaultdict( list )
590      for x in inx.iref_by_sect[i.refid].a['requestVar']:
591           if inx.uid[x].priority <= pmax:
592              e[i.uid].add( inx.uid[x].vid )
593
594              if verbose:
595                cmv = inx.uid[inx.uid[x].vid]
596                if cmv._h.label == 'CMORvar':
597                  si[ cmv.mipTable ].append( inx.uid[x].label )
598#
599# for each variable, calculate the maximum number of years across all the request links which reference that variable.
600##
601## for each request item we have nymax, nenmax, nexmax.
602##
603    nymg = collections.defaultdict( dict )
604##
605## if dataset count rather than volume is wanted, use item 3 from rqiExp tuple.
606    if adsCount:
607      irqi = 3
608    else:
609      irqi = 2
610
611    sgg = set()
612    for v in vars:
613      s = set()
614      sg = collections.defaultdict( set )
615      cc2 = collections.defaultdict( set )
616      cc2s = collections.defaultdict( c1s )
617      for i in l1p:
618##################
619        if (exset == None or i.esid in exset) and v in e[i.rlid]:
620          ix = inx.uid[i.esid]
621          rl = inx.uid[i.rlid]
622          sgg.add( rl.grid )
623          if rl.grid in ['100km','1deg','2deg']:
624            grd = rl.grid
625          else:
626            grd = 'native'
627
628          this = None
629          if exset == None:
630            thisz = 100
631##
632## for a single experiment, look up n years, and n ensemble.
633## should have nstart????
634##
635          elif exi._h.label == 'experiment' or ix._h.label == 'experiment':
636            this = None
637            if ex in self.rqiExp[i.uid][1]:
638              this = self.rqiExp[i.uid][1][ex]
639            elif ix.uid in self.rqiExp[i.uid][1]:
640              this = self.rqiExp[i.uid][1][ix.uid]
641            if this != None:
642              thisns = this[-3]
643              thisny = this[-2]
644              thisne = this[-1]
645              cc2s[grd].a[u].add( filter1( thisns*thisny*thisne, i.nymax) )
646          else:
647            thisz = None
648            if 'experiment' in inx.iref_by_sect[i.esid].a:
649              for u in inx.iref_by_sect[i.esid].a['experiment']:
650                if u in self.rqiExp[i.uid][1]:
651                  this = self.rqiExp[i.uid][1][u]
652                  thisns = this[-3]
653                  thisny = this[-2]
654                  thisne = this[-1]
655##
656###   aggregate year count for each experiment and output grid
657## clarify definition and usage of nymax -- should be redundant ... could be replaced by inward references from "timeSlice"
658                  cc2s[grd].a[u].add( filter1( thisns*thisny*thisne, i.nymax) )
659
660          if exset != None:
661            sg[grd].add( self.rqiExp[i.uid][irqi] )
662     
663###
664### sum over experiments of maximum within each experiment
665###
666      for g in sg:
667        nymg[v][g] = sum( [max( cc2s[g].a[k] ) for k in cc2s[g].a] )
668
669    szv = {}
670    ov = []
671    for v in vars:
672      if 'requestVar' not in inx.iref_by_sect[v].a:
673         print ( 'Variable with no request ....: %s, %s' % (inx.uid[v].label, inx.uid[v].mipTable) )
674      try:
675        szv[v] = self.sz[inx.uid[v].stid]*npy[inx.uid[v].frequency]
676      except:
677        if inx.uid[v].stid not in self.sz:
678          print ('ERROR: size not found for stid %s (v=%s, %s)' % (inx.uid[v].stid,v,inx.uid[v].label) )
679        if inx.uid[v].frequency not in npy:
680          print ('ERROR: npy not found for frequency %s (v=%s, %s)' % (inx.uid[v].frequency,v,inx.uid[v].label) )
681        szv[v] = 0
682      ov.append( self.dq.inx.uid[v] )
683
684    ff = {}
685    for v in vars:
686      if adsCount:
687        ff[v] = 1
688      else:
689        if 'native' in nymg[v]:
690          ff[v] = szv[v]
691          ny = nymg[v]['native']
692        else:
693          if len( nymg[v] ) > 1:
694            print ( '########### Selecting first in list .............' )
695          ks0 = nymg[v].keys()
696          if len(ks0) == 0:
697            ff[v] = 0.
698            ny = 0.
699          else:
700            ks = list( nymg[v].keys() )[0]
701            ny = nymg[v][ks]
702            if inx.uid[v].stid in self.szg[ks]:
703              ff[v] = self.szg[ks][ inx.uid[v].stid ] * npy[inx.uid[v].frequency]
704            else:
705              ff[v] = 0.
706
707        if inx.uid[v].frequency != 'monClim':
708          ff[v] = ff[v]*ny
709
710    ee = self.listIndexDual( ov, 'mipTable', 'label', acount=None, alist=None, cdict=ff, cc=cc )
711    self.ngptot = sum( [  ff[v]  for v in vars] )
712    return (self.ngptot, ee, ff )
713
714  def esid_to_exptList(self,esid,deref=False,full=False):
715    if not esid in self.dq.inx.uid:
716      mlg.prnt ( 'Attempt to dereferece invalid uid: %s' % esid )
717      raise
718
719    if self.dq.inx.uid[esid]._h.label == 'experiment':
720      expts = [esid,]
721    elif self.dq.inx.uid[esid]._h.label != 'remarks':
722      if esid in self.dq.inx.iref_by_sect and 'experiment' in self.dq.inx.iref_by_sect[esid].a:
723        expts = list( self.dq.inx.iref_by_sect[esid].a['experiment'][:] )
724      else:
725        expts = []
726
727## add in groups and mips for completeness
728##
729      if full:
730        if self.dq.inx.uid[esid]._h.label == 'mip':
731          s = set()
732          for e in expts:
733            if self.dq.inx.uid[e]._h.label != 'experiment':
734              mlg.prnt ( 'ERROR: %s, %s, %s ' % (esid,e, self.dq.inx.uid[e].title ) )
735            s.add( self.dq.inx.uid[e].egid )
736          for i in s:
737            expts.append( i )
738        expts.append( esid )
739    else:
740      return None
741
742    if self.tierMax > 0:
743      expts1 = []
744      for i in expts:
745        if self.dq.inx.uid[i]._h.label == 'experiment':
746          if self.dq.inx.uid[i].tier[0] <= self.tierMax:
747            expts1.append( i )
748        elif self.dq.inx.uid[i]._h.label == 'exptgroup':
749          if self.dq.inx.uid[i].tierMin <= self.tierMax:
750            expts1.append( i )
751        else:
752            expts1.append( i )
753    else:
754      expts1 = expts
755
756    if deref:
757      return [self.dq.inx.uid[e] for e in expts1]
758    else:
759      return expts1
760##
761## need to call this on load
762## then use instead of i.ny etc below
763##
764  def requestItemExpAll( self ):
765    self.rqiExp = {}
766    for rqi in self.dq.coll['requestItem'].items:
767      a,b,c,d,e = self.requestItemExp( rqi )
768      if a != None:
769        self.rqiExp[rqi.uid] = (a,b,c,d,e)
770
771  def requestItemExp( self, rqi ):
772    assert rqi._h.label == "requestItem", 'Argument to requestItemExp must be a requestItem'
773    tsl = None
774    if 'tslice' in rqi.__dict__:
775      ts = self.dq.inx.uid[ rqi.tslice ]
776      if ts._h.label == 'timeSlice':
777        if ts.type == 'simpleRange':
778          tsl = (ts.label,'simpleRange', ts.start,ts.end)
779        elif ts.type == 'branchedYears':
780          tsl = (ts.label,'%s:%s' % (ts.type,ts.child), ts.start,ts.end)
781        else:
782          tsl = (ts.label, ts.type, None, None )
783
784     
785    u = rqi.esid
786    if self.dq.inx.uid[u]._h.label == 'experiment':
787      expts = [u,]
788    elif self.dq.inx.uid[u]._h.label != 'remarks':
789      if u in self.dq.inx.iref_by_sect and 'experiment' in self.dq.inx.iref_by_sect[u].a:
790        expts = self.dq.inx.iref_by_sect[u].a['experiment']
791      else:
792        expts = []
793    else:
794      return (None, None, None, None,None)
795
796    if self.tierMax > 0:
797      expts = [i for i in expts if self.dq.inx.uid[i].tier[0] <= self.tierMax]
798
799    self.multiTierOnly = False
800    if self.multiTierOnly:
801      expts = [i for i in expts if len(self.dq.inx.uid[i].tier) > 1]
802      print ('Len expts: %s' % len(expts) )
803
804    if len(expts) > 0:
805      e = [self.dq.inx.uid[i] for i in expts]
806      for i in e:
807        if i._h.label != 'experiment':
808          mlg.prnt ( 'ERROR: %s, %s, %s ' % ( u,i._h.label, i.label, i.title ) )
809      dat2 = {}
810      for i in e:
811        dat2[i.uid] = (i.ntot, i.yps, i.ensz, i.tier, i.nstart, filter1(i.yps,rqi.nymax), filter2(i.ensz,rqi.nenmax,i.tier,self.tierMax) )
812
813      nytot = sum( [dat2[x][-2]*dat2[x][-3] for x in dat2 ] )
814      netot = sum( [dat2[x][-1] for x in dat2 ] )
815    else:
816      dat2 = {}
817      nytot = 0
818      netot = 0
819   
820##
821## to get list of years per expt for each requestLink ... expts is union of all dat2 keys,
822## and want max of dat2[x][0] for each experiment x.
823##
824    return (expts, dat2, nytot, netot, tsl )
825
826  def setTierMax( self, tierMax ):
827    """Set the maxium tier and recompute request sizes"""
828    if tierMax != self.tierMax:
829      self.tierMax = tierMax
830      self.requestItemExpAll(  )
831
832  def summaryByMip( self, pmax=1 ):
833    bytesPerFloat = 2.
834    for m in self.mipls:
835      v = self.volByMip( m, pmax=pmax )
836      mlg.prnt ( '%12.12s: %6.2fTb' % (m,v*bytesPerFloat*1.e-12) )
837
838  def rqlByMip( self, mip):
839    if mip == 'TOTAL':
840        mip = self.mips
841
842    if type(mip) in [type( '' ),type( u'') ]:
843      if mip not in self.mips:
844        mlg.prnt ( self.mips )
845        raise baseException( 'rqiByMip: Name of mip not recognised: %s' % mip )
846      l1 = [i for i in  self.dq.coll['requestLink'].items if i.mip == mip]
847    elif type(mip) in [ type( set()), type( [] ) ]:
848      nf = [ m for m in mip if m not in self.mips]
849      if len(nf) > 0:
850          raise baseException( 'rqlByMip: Name of mip(s) not recognised: %s' % str(nf) )
851      l1 = [i for i in  self.dq.coll['requestLink'].items if i.mip in mip]
852    elif type(mip) == type( dict()):
853      nf = [ m for m in mip if m not in self.mips]
854      if len(nf) > 0:
855        raise baseException( 'rqlByMip: Name of mip(s) not recognised: %s' % str(nf) )
856      l1 = []
857      for i in  self.dq.coll['requestLink'].items:
858        if i.mip in mip:
859          ok = False
860          if len( mip[i.mip] ) == 0:
861            ok = True
862          else:
863            for ol in self.dq.inx.iref_by_sect[i.uid].a['objectiveLink']:
864              o = self.dq.inx.uid[ol]
865              if self.dq.inx.uid[o.oid].label in mip[i.mip]:
866                ok = True
867          if ok:
868              l1.append( i )
869    else:
870      raise baseException( 'rqiByMip: "mip" (1st explicit argument) should be type string or set: %s -- %s' % (mip, type(mip))   )
871
872    return l1
873
874  def rqiByMip( self, mip):
875    l1 = self.rqlByMip( mip )
876    if len(l1) == 0:
877       return []
878    l2 = [] 
879    for i in l1:
880       if 'requestItem' in self.dq.inx.iref_by_sect[i.uid].a:
881          for u in self.dq.inx.iref_by_sect[i.uid].a['requestItem']:
882               l2.append( self.dq.inx.uid[u] )
883
884    l20 = self.rqiByMip0( mip )
885    for i in l20:
886      assert i in l2
887    return l2
888   
889   
890  def rqiByMip0( self, mip):
891
892    if mip == 'TOTAL':
893        mip = self.mips
894    if type(mip) in [type( '' ),type( u'') ]:
895      if mip not in self.mips:
896        mlg.prnt ( self.mips )
897        raise baseException( 'rqiByMip: Name of mip not recognised: %s' % mip )
898      l1 = [i for i in  self.dq.coll['requestItem'].items if i.mip == mip]
899    elif type(mip) in [ type( set()), type( [] ) ]:
900      nf = [ m for m in mip if m not in self.mips]
901      if len(nf) > 0:
902          raise baseException( 'rqiByMip: Name of mip(s) not recognised: %s' % str(nf) )
903      l1 = [i for i in  self.dq.coll['requestItem'].items if i.mip in mip]
904    elif type(mip) == type( dict()):
905      nf = [ m for m in mip if m not in self.mips]
906      if len(nf) > 0:
907        raise baseException( 'rqiByMip: Name of mip(s) not recognised: %s' % str(nf) )
908      l1 = []
909      for i in  self.dq.coll['requestLink'].items:
910        if i.mip in mip:
911          ok = False
912          if len( mip[i.mip] ) == 0:
913            ok = True
914          else:
915            for ol in self.dq.inx.iref_by_sect[i.uid].a['objectiveLink']:
916              o = self.dq.inx.uid[ol]
917              if self.dq.inx.uid[o.oid].label in mip[i.mip]:
918                ok = True
919          if ok:
920              if 'requestItem' in self.dq.inx.iref_by_sect[i.uid].a:
921                for u in self.dq.inx.iref_by_sect[i.uid].a['requestItem']:
922                  l1.append( self.dq.inx.uid[u] )
923    else:
924      raise baseException( 'rqiByMip: "mip" (1st explicit argument) should be type string or set: %s -- %s' % (mip, type(mip))   )
925
926    return l1
927
928  def checkDir(self,odir,msg):
929      if not os.path.isdir( odir ):
930         try:
931            os.mkdir( odir )
932         except:
933            print ('\n\nFailed to make directory "%s" for: %s: make necessary subdirectories or run where you have write access' % (odir,msg) )
934            print ( '\n\n' )
935            raise
936         print ('Created directory %s for: %s' % (odir,msg) )
937
938  def xlsByMipExpt(self,m,ex,pmax,odir='xls',xls=True,txt=False,txtOpts=None):
939    import scope_utils
940    mxls = scope_utils.xlsTabs(self,tiermax=self.tierMax,pmax=pmax,xls=xls, txt=txt, txtOpts=txtOpts,odir=odir)
941    mlab = makeTables.setMlab( m )
942    mxls.run( m, mlab=mlab )
943
944  def cmvByInvMip( self, mip,pmax=1,includeYears=False, exptFilter=None,exptFilterBlack=None ):
945    mips = set( self.mips[:] )
946    if type(mip) == type( '' ):
947        mips.discard( mip )
948    else:
949      for m in mip:
950        mips.discard( m )
951
952    return self.cmvByMip( mips,pmax=pmax,includeYears=includeYears, exptFilter=exptFilter, exptFilterBlack=exptFilterBlack )
953
954  def cmvByMip( self, mip,pmax=1,includeYears=False, exptFilter=None, exptFilterBlack=None ):
955    if exptFilter != None:
956      assert type(exptFilter) == type( set() ), 'Argument exptFilter must be None or a set: %s' % str(exptFilter)
957    if exptFilterBlack != None:
958      assert type(exptFilterBlack) == type( set() ), 'Argument exptFilterBlack must be None or a set: %s' % str(exptFilterBlack)
959      if exptFilter != None:
960        assert len( exptFilter.difference( exptFilterBlack ) ) > 0, 'If exptFilter and exptFilterBlack are both set, exptFilter must have non-black listed elements' 
961
962    l1,ee = self.rvgByMip( mip, includePreset=True, returnLinks=True )
963    if includeYears:
964      expys = self.exptYears( l1, ex=exptFilter, exBlack=exptFilterBlack )
965      cc = collections.defaultdict( set )
966      ccts = collections.defaultdict( set )
967    ss = set()
968    for pr in ee:
969### loop over request  var groups.
970      for i in ee[pr]:
971        if 'requestVar' in self.dq.inx.iref_by_sect[i.uid].a:
972#
973# loop over request vars in group
974#
975          for x in self.dq.inx.iref_by_sect[i.uid].a['requestVar']:
976            i1 = self.dq.inx.uid[x]
977
978            if (pr == -1 and i1.priority <= pmax) or (pr > 0 and pr <= pmax):
979              if includeYears and i1.vid in self.cmvGridId:
980                ##assert i.uid in expys, 'No experiment info found for requestVarGroup: %s' % i.uid
981                ## may have no entry as a consequence of tierMin being set in the requestLink(s).
982                assert i1.vid in self.cmvGridId, 'No grid identification lookup found for %s: %s' % (i1.label,i1.vid)
983                assert self.cmvGridId[i1.vid] in ['a','o','si','li'], 'Unexpected grid id: %s: %s:: %s' % (i1.label,i1.vid, self.cmvGridId[i1.vid])
984                gflg = {'si':'','li':''}.get( self.cmvGridId[i1.vid], self.cmvGridId[i1.vid] )
985                rtl = True
986
987                if i.uid in expys:
988                  if rtl:
989                    for e,grd in expys[i.uid]:
990                        if exptFilter == None or e in exptFilter:
991                          if grd == 'DEF':
992                            if gflg == 'o':
993                              grd1 = '1deg'
994                            else:
995                              grd1 = 'native'
996                          else:
997                            grd1 = grd
998                          cc[(i1.vid,e,grd1)].add( expys[i.uid][e,grd] )
999                          if i.uid in self.tsliceDict and e in self.tsliceDict[i.uid]:
1000                            ccts[(i1.vid,e)].add( self.tsliceDict[i.uid][e] )
1001
1002                  else:
1003
1004                   for gf in expys[i.uid]:
1005                    for e,grd in expys[i.uid][gf]:
1006                      if grd in ["1deg",'2deg'] or gf == gflg:
1007                        if exptFilter == None or e in exptFilter:
1008                          cc[(i1.vid,e,grd)].add( expys[i.uid][gf][e,grd] )
1009              else:
1010                print ( 'SKIPPING %s: %s' % (i1.label,i1.vid) )
1011                ss.add( i1.vid )
1012    if includeYears:
1013      l2 = collections.defaultdict( dict )
1014      l2x = collections.defaultdict( dict )
1015##
1016## this removes lower ranked grids .... but for some groups want different grids for different variable categories
1017##
1018      if self.gridPolicyTopOnly:
1019        for v,e,g in cc:
1020          l2x[(v,e)][g] = max( list( cc[(v,e,g)] ) )
1021        for v,e in l2x:
1022          if len( l2x[(v,e)].keys() ) == 1:
1023             g,val = list( l2x[(v,e)].items() )[0]
1024          else:
1025            if 'native' in l2x[(v,e)].keys():
1026               g = 'native'
1027               val = l2x[(v,e)][g]
1028            else:
1029               g = sorted( list( l2x[(v,e)].keys() ) )[0]
1030               val = l2x[(v,e)][g]
1031          l2[v][(e,g)] = val
1032      else:
1033        for v,e,g in cc:
1034          l2[v][(e,g)] = max( list( cc[(v,e,g)] ) )
1035
1036      l2ts = collections.defaultdict( dict )
1037      for v in l2:
1038        for e,g in l2[v]:
1039          if (v,e) in ccts:
1040            if len( ccts[(v,e)] ) > 1:
1041              rc, ts, msg = sortTimeSlice( ccts[(v,e)] )
1042              if rc == 1:
1043                l2ts[v][e] = ts
1044              elif rc == 2:
1045                yl = range( ts[0][2], ts[0][3] + 1) + range( ts[1][2], ts[1][3] + 1)
1046                l2ts[v][e] = ('_union', 'YEARLIST', len(yl), yl )
1047              else:
1048                print ('TIME SLICE MULTIPLE OPTIONS FOR : %s, %s, %s, %s' % (v,e,str(ccts[(v,e)]), msg ) )
1049            else:
1050              l2ts[v][e] = ccts[(v,e)].pop()
1051      return l2, l2ts
1052    else:
1053      l2 = sorted( [i for i in [self.dq.inx.uid[i] for i in ss] if i._h.label != 'remarks'], key=lambda x: x.label )
1054      return l2
1055
1056  def exptFilterList(self,val,option,ret='uid'):
1057    if type( val ) not in [[],()]:
1058      val = [val,]
1059
1060    if option == 'lab':
1061      v0 = val[:]
1062      val = []
1063      mm = []
1064      for v in v0:
1065        if v not in self.exptByLabel:
1066          mm.append( v )
1067        else:
1068          val.append( self.exptByLabel[v] )
1069
1070      assert len(mm) == 0, 'Experiment names not all recognised: %s' % str(mm)
1071
1072    oo = set()
1073    for v in val:
1074      i = self.dq.inx.uid[v]
1075      if i._h.label in ['exptgroup','mip']:
1076        if 'experiment' in self.dq.inx.iref_by_sect[i.uid].a:
1077          for u in self.dq.inx.iref_by_sect[i.uid].a['experiment']:
1078            oo.add( u )
1079      elif i._h.label == 'experiment':
1080            oo.add( i.uid )
1081      else:
1082        print ('WARNING .. skipping request for experiment which links to record of type %s' % i._h.label )
1083    return oo
1084   
1085  def getFreqStrSummary(self,mip,pmax=1):
1086##
1087## get a dictionary keyed on CMORvar uid, containing dictionary keyed on (experiment, grid) with value as number of years.
1088##
1089    if not self.uniqueRequest:
1090      cmv, self.cmvts = self.cmvByMip(mip,pmax=pmax,includeYears=True,exptFilter=self.exptFilter,exptFilterBlack=self.exptFilterBlack)
1091    else:
1092      cmv1, cmvts1 = self.cmvByInvMip(mip,pmax=pmax,includeYears=True,exptFilter=self.exptFilter,exptFilterBlack=self.exptFilterBlack)
1093      cmv2, cmvts2 = self.cmvByMip('TOTAL',pmax=pmax,includeYears=True,exptFilter=self.exptFilter,exptFilterBlack=self.exptFilterBlack)
1094      cmv = self.differenceSelectedCmvDict(  cmv1, cmv2 )
1095 
1096    self.selectedCmv = cmv
1097    return self.cmvByFreqStr( cmv )
1098
1099  def differenceSelectedCmvDict( self, cmv1, cmv2 ):
1100      """Return the diffence between two dictionaries of cmor variables returned by self.cmvByMip.
1101         The dictionaries contain dictionaries of values. Differences may be subdictionaries not present,
1102         elements of sub-dictionaries not present, or elements of sub-dictionaries present with different values.
1103         A one sided difference is returned."""
1104
1105      cmv = {}
1106      for i in cmv2:
1107        if i not in cmv1:
1108          cmv[i] = cmv2[i]
1109        else:
1110          eei = {}
1111          for t in cmv2[i]:
1112            if t not in cmv1[i]:
1113              eei[t] = cmv2[i][t]
1114            else:
1115              if cmv2[i][t] > cmv1[i][t]:
1116                 eei[t] = cmv2[i][t] - cmv1[i][t]
1117          if len( eei.keys() ) != 0:
1118            cmv[i] = eei
1119      return cmv
1120
1121  def cmvByFreqStr(self,cmv,asDict=True,exptFilter=None,exptFilterBlack=None):
1122    if exptFilter != None:
1123      assert type(exptFilter) == type( set() ), 'Argument exptFilter must be None or a set: %s' % str(exptFilter)
1124    if exptFilterBlack != None:
1125      assert type(exptFilterBlack) == type( set() ), 'Argument exptFilterBlack must be None or a set: %s' % str(exptFilterBlack)
1126      if exptFilter != None:
1127        assert len( exptFilter.difference( exptFilterBlack ) ) > 0, 'If exptFilter and exptFilterBlack are both set, exptFilter must have non-black listed elements' 
1128
1129    cc = collections.defaultdict( list )
1130    for i in cmv:
1131      if asDict:
1132        ii = self.dq.inx.uid[i]
1133        if ii._h.label != 'remarks':
1134          st = self.dq.inx.uid[ ii.stid ]
1135          cc0 = collections.defaultdict( float )
1136          cc1 = collections.defaultdict( int )
1137          se = collections.defaultdict( set )
1138          for e,g in cmv[i]:
1139            cc0[g] += cmv[i][(e,g)]
1140            cc1[g] += 1
1141            se[g].add(e)
1142          for g in cc0:
1143            g1 = 'native'
1144            if self.isLatLon[st.spid] != False:
1145              g1 = g
1146              if g1 == 'DEF' and self.isLatLon[st.spid] == 'o':
1147                  g1 = '1deg'
1148              else:
1149                  g1 = 'native'
1150            g1 = g
1151
1152            cc[ (st.spid,st.odims,ii.frequency,g1) ].append( (i,cc0[g],cc1[g],se[g]) )
1153
1154      else:
1155        st = self.dq.inx.uid[ i.stid ]
1156        cc[ (st.spid,st.odims,i.frequency) ].append( i.label )
1157
1158    self.thiscmvset = set()
1159    c2 = collections.defaultdict( dict )
1160    sf = set()
1161    if asDict:
1162      for s,o,f,g in cc.keys():
1163        c2[(s,o,g)][f] = cc[ (s,o,f,g) ]
1164        sf.add( f )
1165    else:
1166      for s,o,f in cc.keys():
1167        c2[(s,o)][f] = cc[ (s,o,f) ]
1168        sf.add( f )
1169    lf = sorted( list(sf) )
1170    c3 = collections.defaultdict( dict )
1171
1172    for tt in sorted( c2.keys() ):
1173      if asDict:
1174        s,o,g = tt
1175      else:
1176        s,o = tt
1177        g = 'native'
1178      i = self.dq.inx.uid[ s ]
1179
1180      if asDict:
1181        for f in c2[tt]:
1182            isClim = f.lower().find( 'clim' ) != -1
1183            ny = 0
1184            expts = set()
1185            labs = []
1186            labs = collections.defaultdict( int )
1187            ccx = collections.defaultdict( list )
1188            for cmvi, ny1, ne, eset in c2[tt][f]:
1189              ccx[cmvi].append( (ny1, ne, eset) )
1190            net = 0
1191            for cmvi in ccx:
1192              if len( ccx[cmvi] ) == 1:
1193                 ny1, ne, eset = ccx[cmvi][0]
1194              else:
1195                 ny1, ne, eset = ( 0,0,set() )
1196                 for a,b,s in ccx[cmvi]:
1197                   ny1 += a
1198                   ne += b
1199                   eset = eset.union(  s )
1200             
1201              net += ne
1202              if len(eset) != ne:
1203                print ( 'WARNING: inconsistency in volume estimate ... possible duplication for %s,%s' % (cmvi,f) )
1204              for e in eset:
1205                elab = self.dq.inx.uid[e].label
1206                expts.add(elab)
1207
1208              if exptFilter != None:
1209                expts = exptFilter.intersection( expts )
1210              if exptFilterBlack != None:
1211                expts = expts.difference( exptFilterBlack )
1212
1213              if len(expts) > 0:
1214                lab = self.dq.inx.uid[cmvi].label
1215                self.thiscmvset.add( cmvi )
1216                ny += ny1
1217                labs[cmvi] += ny1
1218            ne = len( expts )
1219            nn = len( labs.keys() )
1220             
1221            if isClim:
1222              ny = net/float(nn)
1223            else:
1224              ny = ny/float(nn)
1225            assert tt[2] in ['native','1deg','2deg'], 'BAD grid identifier: %s' % str(tt)
1226            c3[tt][f] = (nn,ny,ne, labs,expts)
1227    return (sf,c3)
1228
1229  def getStrSz( self, g, stid=None, s=None, o=None, tt=False ):
1230    assert stid == None or (s==None and o==None), 'Specify either stid or s and o'
1231    assert stid != None or (s!=None and o!=None), 'Specify either stid or s and o'
1232
1233    if stid != None:
1234      st = self.dq.inx.uid[stid]
1235      if st._h.label != 'remarks':
1236        s = st.spid
1237        o = st.odims
1238      else:
1239        self.strSz[ (stid,g) ] = (False,0)
1240        if tt:
1241          return (self.strSz[ (stid,g) ], None)
1242        else:
1243          return self.strSz[ (stid,g) ]
1244
1245    g1 = g
1246    if g1 == 'DEF':
1247          if self.isLatLon[s] == 'o':
1248             g1 = '1deg'
1249          else:
1250             g1 = 'native'
1251    if (s,o,g) not in self.strSz:
1252
1253        if o == '':
1254           sf = 1
1255        elif o in self.odsz:
1256           sf = self.odsz[o]
1257        else:
1258           sf = 5
1259
1260        if type( sf ) == type( () ):
1261           sf = sf[0]
1262
1263
1264        try:
1265          if g1 != 'native' and self.isLatLon[s] != False:
1266            szg = self.szgss[g1][s]
1267          else:
1268            szg = self.szss[s]
1269        except:
1270          print ( 'Failed to get size for: %s, %s, %s' % (g,g1,s ) )
1271          raise
1272
1273        szg = szg * sf
1274        self.strSz[ (s,o,g) ] = (True,szg)
1275
1276    if tt:
1277      return (self.strSz[ (s,o,g) ], (s,o,g1) )
1278    else:
1279      return self.strSz[ (s,o,g) ]
1280
1281  def rvgByMip( self, mip, years=False, includePreset=False, returnLinks=False ):
1282    l1 = self.rqlByMip( mip )
1283    if includePreset:
1284      cc = collections.defaultdict( set )
1285      ss = set()
1286      for i in l1:
1287        if 'requestItem' in self.dq.inx.iref_by_sect[i.uid].a:
1288          prs = set()
1289          for x in self.dq.inx.iref_by_sect[i.uid].a['requestItem']:
1290             prs.add(self.dq.inx.uid[x].preset)
1291
1292          for p in prs:
1293            assert p in [-1,1,2,3], 'Bad preset value'
1294            cc[p].add( i.refid )
1295      ee = {}
1296      for p in cc:
1297        l2 = sorted( [self.dq.inx.uid[i] for i in cc[p]], key=lambda x: x.label )
1298        ee[p] = l2
1299      if returnLinks:
1300        return (l1,ee)
1301      else:
1302        return ee
1303    else:
1304      ss = set( [i.refid for i in l1] )
1305      l2 = sorted( [self.dq.inx.uid[i] for i in ss], key=lambda x: x.label )
1306      if returnLinks:
1307        return (l1,l2)
1308      else:
1309        return l2
1310
1311  def volByMip( self, mip, pmax=2, retainRedundantRank=False, intersection=False, adsCount=False, exptid=None):
1312
1313    l1 = self.rqiByMip( mip )
1314     
1315    #### The set of experiments/experiment groups:
1316    if exptid == None:
1317      exps = self.experiments
1318    elif type( exptid ) == type(''):
1319      exps = set( [exptid,] )
1320    else:
1321      assert type( exptid ) == type( set() ),'exptid arg to volByMip must be None, string or set: %s' % type( exptid )
1322      exps = exptid
1323   
1324    self.volByE = {}
1325    vtot = 0
1326    cc = collections.defaultdict( col_count )
1327    self.allVars = set()
1328    for e in exps:
1329      expts = self.esid_to_exptList(e,deref=True,full=False)
1330      if expts not in  [None,[]]:
1331        for ei in expts:
1332          self.volByE[ei.label] = self.volByExpt( l1, ei.uid, pmax=pmax, cc=cc, retainRedundantRank=retainRedundantRank, intersection=intersection, adsCount=adsCount )
1333          vtot += self.volByE[ei.label][0]
1334        self.allVars = self.allVars.union( self.vars )
1335    self.indexedVol = cc
1336
1337    return vtot
1338
1339  def listIndexDual(self, ll, a1, a2, acount=None, alist=None, cdict=None, cc=None ):
1340    do_count = acount != None
1341    do_list = alist != None
1342    assert not (do_count and do_list), 'It is an error to request both list and count'
1343    if not (do_count or do_list):
1344      acount = '__number__'
1345      do_count = True
1346
1347    if cc == None:
1348      if do_count:
1349        cc = collections.defaultdict( col_count )
1350      elif do_list:
1351        cc = collections.defaultdict( col_list )
1352
1353    if do_count:
1354      for l in ll:
1355        if cdict != None:
1356          v = cdict[l.uid]
1357        elif acount == '__number__':
1358          v = 1
1359        else:
1360          v = l.__dict__[acount]
1361
1362        cc[ l.__dict__[a1] ].a[ l.__dict__[a2] ] += v
1363    elif do_list:
1364      for l in ll:
1365        if cdict != None:
1366          v = cdict[l.uid]
1367        elif alist == '__item__':
1368          v = l
1369        else:
1370          v = l.__dict__[alist]
1371        cc[ l.__dict__[a1] ].a[ l.__dict__[a2] ].append( v )
1372
1373    od = {}
1374    for k in cc.keys():
1375      d2 = {}
1376      for k2 in cc[k].a.keys():
1377        d2[k2] = cc[k].a[k2]
1378      od[k] = d2
1379    return od
1380
1381class dreqUI(object):
1382  """Data Request Command line.
1383-------------------------
1384      -v : print version and exit;
1385      --unitTest : run some simple tests;
1386      -m <mip>:  MIP of list of MIPs (comma separated; for objective selection see note [1] below);
1387      -l <options>: List for options:
1388              o: objectives
1389              e: experiments
1390      -q <options>: List information about the schema:
1391              s: sections
1392              <section>: attributes for a section
1393              <section:attribute>: definition of an attribute.
1394      -h :       help: print help text;
1395      -e <expt>: experiment;
1396      -t <tier> maxmum tier;
1397      -p <priority>  maximum priority;
1398      --xls : Create Excel file with requested variables;
1399      --sf : Print summary of variable count by structure and frequency;
1400      --xfr : Output variable lists in sheets organised by frequency and realm instead of by MIP table;
1401      --SF : Print summary of variable count by structure and frequency for all MIPs;
1402      --grdpol <native|1deg> :  policy for default grid, if MIPs have not expressed a preference;
1403      --allgrd :  When a variable is requested on multiple grids, archive all grids requested (default: only the finest resolution);
1404      --unique :  List only variables which are requested uniquely by this MIP, for at least one experiment;
1405      --txt : Create text file with requested variables;
1406      --mcfg : Model configuration: 7 integers, comma separated, 'nho','nlo','nha','nla','nlas','nls','nh1'
1407                 default: 259200,60,64800,40,20,5,100
1408      --txtOpts : options for content of text file: (v|c)[(+|-)att1[,att2[...]]]
1409      --xlsDir <directory> : Directory in which to place variable listing [xls];
1410      --printLinesMax <n>  : Maximum number of lines to be printed (default 20)
1411      --printVars    : If present, a summary of the variables (see --printLinesMax) fitting the selection options will be printed
1412      --intersection : Analyse the intersection of requests rather than union.
1413
1414NOTES
1415-----
1416[1] A set of objectives within a MIP can be specified in the command line. The extended syntax of the "-m" argument is:
1417-m <mip>[:objective[.obj2[.obj3 ...]]][,<mip2]...]
1418
1419e.g.
1420drq -m HighResMIP:Ocean.DiurnalCycle
1421"""
1422  def __init__(self,args):
1423    self.adict = {}
1424    self.knownargs = {'-m':('m',True), '-p':('p',True), '-e':('e',True), '-t':('t',True), \
1425                      '-h':('h',False), '--printLinesMax':('plm',True), \
1426                      '-l':('l',True),
1427                      '-q':('q',True),
1428                      '--printVars':('vars',False), '--intersection':('intersection',False), \
1429                      '--count':('count',False), \
1430                      '--txt':('txt',False), \
1431                      '--sf':('sf',False), \
1432                      '--xfr':('xfr',False), \
1433                      '--SF':('SF',False), \
1434                      '--grdpol':('grdpol',True), \
1435                      '--allgrd':('allgrd',False), \
1436                      '--unique':('unique',False), \
1437                      '--mcfg':('mcfg',True), \
1438                      '--txtOpts':('txtOpts',True), \
1439                      '--xlsDir':('xlsdir',True), '--xls':('xls',False) \
1440                       } 
1441    aa = args[:]
1442    notKnownArgs = []
1443    while len(aa) > 0:
1444      a = aa.pop(0)
1445      if a in self.knownargs:
1446        b = self.knownargs[a][0]
1447        if self.knownargs[a][1]:
1448          v = aa.pop(0)
1449          self.adict[b] = v
1450        else:
1451          self.adict[b] = True
1452      else:
1453        notKnownArgs.append(a)
1454
1455    assert self.checkArgs( notKnownArgs ), 'FATAL ERROR 001: Arguments not recognised: %s' % (str(notKnownArgs) )
1456
1457    if 'm' in self.adict:
1458      if self.adict['m'] == '_all_':
1459        pass
1460      elif self.adict['m'].find( ':' ) != -1:
1461        ee = {}
1462        for i in self.adict['m'].split(','):
1463          bits =  i.split( ':' )
1464          if len( bits ) == 1:
1465             ee[bits[0]] = []
1466          else:
1467             assert len(bits) == 2, 'Cannot parse %s' % self.adict['m']
1468             ee[bits[0]] = bits[1].split( '.' )
1469        self.adict['m'] = ee
1470      else:
1471        self.adict['m'] = set(self.adict['m'].split(',') )
1472
1473    if 'grdpol' in self.adict:
1474      assert self.adict['grdpol'] in ['native','1deg'], 'Grid policy argument --grdpol must be native or 1deg : %s' % self.adict['grdpol']
1475
1476    integerArgs = set( ['p','t','plm'] )
1477    for i in integerArgs.intersection( self.adict ):
1478      self.adict[i] = int( self.adict[i] )
1479
1480    self.intersection = self.adict.get( 'intersection', False )
1481
1482 
1483  def checkArgs( self, notKnownArgs ):
1484    if len( notKnownArgs ) == 0:
1485      return True
1486    print ('--------------------------------------')
1487    print ('------------  %s Arguments Not Recognised ------------' % len(notKnownArgs) )
1488    k = 0
1489    for x in notKnownArgs:
1490      k += 1
1491      if x[1:] in self.knownargs:
1492        print ( '%s PERHAPS %s instead of %s' % (k, x[1:],x) )
1493      elif '-%s' % x in self.knownargs:
1494        print ( '%s PERHAPS -%s instead of %s' % (k, x,x) )
1495      elif x[0] == '\xe2':
1496        print ( '%s POSSIBLY -- (double hyphen) instead of long dash in %s' % (k, x) )
1497    print ('--------------------------------------')
1498
1499    return len( notKnownArgs ) == 0
1500     
1501  def run(self, dq=None):
1502    if 'h' in self.adict:
1503      mlg.prnt ( self.__doc__ )
1504      return
1505
1506    if 'q' in self.adict:
1507      if dq == None:
1508        dq = dreq.loadDreq(configOnly=True)
1509      s = self.adict['q']
1510      if self.adict['q'] == 's':
1511        ss = sorted( [(i.title,i.label) for i in dq.coll['__sect__'].items] )
1512        for s in ss:
1513          mlg.prnt( '%16s:: %s' % (s[1],s[0]) )
1514      else:
1515        ss = [i.label for i in dq.coll['__sect__'].items]
1516        if s.find( ':' ) != -1:
1517          s,a = s.split( ':' )
1518        else:
1519          a = None
1520        if s not in ss:
1521          mlg.prnt( 'ERROR: option must be a section; use "-q s" to list sections' )
1522        elif a == None:
1523          x = [i for i in dq.coll['__sect__'].items if i.label == s]
1524          s1 = [i for i in  dq.coll['__main__'].items if 'ATTRIBUTE::%s' % s in i.uid]
1525          mlg.prnt( x[0].title )
1526          mlg.prnt( ' '.join( sorted  ([i.label for i in s1] ) ))
1527        else:
1528          x = [i for i in dq.coll['__main__'].items if i.uid == 'ATTRIBUTE::%s.%s' % (s,a) ]
1529          if len(x) == 0:
1530            mlg.prnt( 'ERROR: attribute not found' )
1531            s1 = [i for i in  dq.coll['__main__'].items if 'ATTRIBUTE::%s' % s in i.uid]
1532            mlg.prnt( 'ATTRIBUTES: ' + ' '.join( sorted  ([i.label for i in s1] ) ))
1533          else:
1534            mlg.prnt( 'Section %s, attribute %s' % (s,a) )
1535            mlg.prnt( x[0].title )
1536            mlg.prnt( x[0].description )
1537      return
1538
1539    if not ('m' in self.adict or 'SF' in self.adict):
1540      mlg.prnt ( 'Current version requires -m or --SF argument'  )
1541      mlg.prnt ( self.__doc__ )
1542      sys.exit(0)
1543
1544    if dq == None:
1545      self.dq = dreq.loadDreq()
1546    else:
1547      self.dq = dq
1548
1549    if 'l' in self.adict:
1550      self.printList()
1551      return
1552
1553    if 'mcfg' in self.adict:
1554      ll = string.split( self.adict['mcfg'], ',' )
1555      assert len(ll) == 7, 'Length of model configuration argument must be 7 comma separated integers: %s' %  self.adict['mcfg']
1556      lli = [ int(x) for x in ll]
1557
1558    self.sc = dreqQuery( dq=self.dq )
1559
1560    if 'grdpol' in self.adict:
1561      self.sc.gridPolicyDefaultNative = self.adict['grdpol'] == 'native'
1562      print ( 'SETTING grid policy: %s' % self.sc.gridPolicyDefaultNative )
1563    if 'allgrd' in self.adict:
1564      self.sc.gridPolicyTopOnly = False
1565      print ( 'SETTING grid policy for multiple preferred grids: %s' % self.sc.gridPolicyTopOnly )
1566    if 'unique' in self.adict:
1567      self.sc.uniqueRequest = True
1568
1569    if 'mcfg' in self.adict:
1570      self.sc.setMcfg( lli )
1571
1572    tierMax = self.adict.get( 't', 1 )
1573    self.sc.setTierMax(  tierMax )
1574    pmax = self.adict.get( 'p', 1 )
1575
1576    makeXls = self.adict.get( 'xls', False )
1577    makeTxt = self.adict.get( 'txt', False )
1578    doSf = 'SF' in self.adict or 'sf' in self.adict
1579    if makeXls or makeTxt or doSf:
1580      xlsOdir = self.adict.get( 'xlsdir', 'xls' )
1581      self.sc.checkDir( xlsOdir, 'xls files' )
1582
1583    tabByFreqRealm = self.adict.get( 'xfr', False )
1584    if 'SF' in self.adict:
1585      self.sc.gridPolicyDefaultNative = True
1586      vs = volsum.vsum( self.sc, odsz, npy, makeTables.makeTab, makeTables.tables, odir=xlsOdir, tabByFreqRealm=tabByFreqRealm )
1587      vs.analAll(pmax)
1588
1589      self.sc.gridPolicyDefaultNative = False
1590      vs = volsum.vsum( self.sc, odsz, npy, makeTables.makeTab, makeTables.tables, odir=xlsOdir, tabByFreqRealm=tabByFreqRealm )
1591      vs.analAll(pmax)
1592
1593      self.sc.setTierMax( 3 )
1594      vs = volsum.vsum( self.sc, odsz, npy, makeTables.makeTab, makeTables.tables, odir=xlsOdir, tabByFreqRealm=tabByFreqRealm )
1595      vs.analAll(3)
1596      return
1597
1598    ok = True
1599    if self.adict['m'] == '_all_':
1600      self.adict['m'] = set(self.sc.mips )
1601      mlab = 'TOTAL'
1602    else:
1603      for i in self.adict['m']:
1604        if i not in self.sc.mips:
1605          ok = False
1606          mlg.prnt ( 'NOT FOUND: %s' % i )
1607      mlab = makeTables.setMlab( self.adict['m'] )
1608    assert ok,'Available MIPs: %s' % str(self.sc.mips)
1609
1610    eid = None
1611    ex = None
1612    if 'e' in self.adict:
1613      ex = self.adict['e']
1614      if ex in self.sc.mipsp:
1615        eid = set( self.dq.inx.iref_by_sect[ex].a['experiment'] )
1616        self.sc.exptFilter = eid
1617      else:
1618        for i in self.dq.coll['experiment'].items:
1619          if i.label == self.adict['e']:
1620            eid = i.uid
1621        assert eid != None, 'Experiment/MIP %s not found' % self.adict['e']
1622        self.sc.exptFilter = set( [eid,] )
1623
1624    ss = set()
1625    for e in ['esm-hist','esm-hist-ext','esm-piControl','piControl-spinup','esm-piControl-spinup']:
1626      ss.add( self.sc.exptByLabel[ e ] )
1627    self.sc.exptFilterBlack = ss
1628
1629    if 'sf' in self.adict:
1630      vs = volsum.vsum( self.sc, odsz, npy, makeTables.makeTab, makeTables.tables, odir=xlsOdir, tabByFreqRealm=tabByFreqRealm )
1631      vs.run( self.adict['m'], 'requestVol_%s_%s_%s' % (mlab,tierMax,pmax), pmax=pmax ) 
1632      vs.anal(olab=mlab,doUnique=False)
1633      mips = self.adict['m']
1634      if type(mips) in [type(set()),type(dict())]:
1635        mips = self.adict['m'].copy()
1636        if len(mips) > 1:
1637          if type(mips) == type(set()):
1638             mips.add( '*TOTAL' )
1639          else:
1640             mips['*TOTAL'] = ''
1641
1642      vs.analAll(pmax,mips=mips,html=False)
1643      ttl = sum( [x for k,x in vs.res['vu'].items()] )*2.*1.e-12
1644      ttl2 = sum( [x for k,x in vs.res['vu'].items()] )*2.*1.e-12
1645      mlg.prnt( 'TOTAL volume: %8.2fTb' % ttl )
1646      return
1647
1648
1649    adsCount = self.adict.get( 'count', False )
1650
1651    self.getVolByMip(pmax,eid,adsCount)
1652    makeXls = self.adict.get( 'xls', False )
1653    makeTxt = self.adict.get( 'txt', False )
1654    if makeXls or makeTxt:
1655      mips = self.adict['m']
1656
1657      if 'txtOpts' in self.adict:
1658        if self.adict['txtOpts'][0] == 'v':
1659          txtOpts = NT_txtopts( 'var' )
1660        else:
1661          txtOpts = NT_txtopts( 'cmv' )
1662      else:
1663        txtOpts=None
1664
1665      self.sc.xlsByMipExpt(mips,eid,pmax,odir=xlsOdir,xls=makeXls,txt=makeTxt,txtOpts=txtOpts)
1666
1667  def printList(self):
1668    mips = self.adict['m']
1669    ee = {}
1670    for i in self.dq.coll['mip'].items:
1671      if i.label in mips:
1672        ee[i.label] = i
1673    if self.adict['l'] in ['o','e']:
1674      targ = {'o':'objective', 'e':'experiment' }[self.adict['l']]
1675      for k in sorted( ee.keys() ):
1676        if targ in self.dq.inx.iref_by_sect[ee[k].uid].a:
1677          for u in self.dq.inx.iref_by_sect[ee[k].uid].a[targ]:
1678            print ( '%s: %s' % (ee[k].label, self.dq.inx.uid[u].label) )
1679    else:
1680      print ('list objective *%s* not recognised (should be e or o)' % self.adict['l'] )
1681     
1682  def getVolByMip(self,pmax,eid,adsCount):
1683
1684    v0 = self.sc.volByMip( self.adict['m'], pmax=pmax, intersection=self.intersection, adsCount=adsCount, exptid=eid )
1685    mlg.prnt ( 'getVolByMip: %s [%s]' % (v0,makeTables.vfmt(v0*2.)) )
1686    cc = collections.defaultdict( int )
1687    for e in self.sc.volByE:
1688      for v in self.sc.volByE[e][2]:
1689          cc[v] += self.sc.volByE[e][2][v]
1690    x = 0
1691    for v in cc:
1692      x += cc[v]
1693   
1694    if python2:
1695      vl = sorted( cc.keys(), cmp=cmpd(cc).cmp, reverse=True )
1696    else:
1697      vl = sorted( cc.keys(), key=lambda x: cc[x], reverse=True )
1698    if self.adict.get( 'vars', False ):
1699      printLinesMax = self.adict.get( 'plm', 20 )
1700      if printLinesMax > 0:
1701        mx = min( [printLinesMax,len(vl)] )
1702      else:
1703        mx = len(vl)
1704
1705      for v in vl[:mx]:
1706        mlg.prnt ( '%s.%s: %s' % (self.dq.inx.uid[v].mipTable,self.dq.inx.uid[v].label, makeTables.vfmt( cc[v]*2. ) ) )
1707      if mx < len(vl):
1708        mlg.prnt ( '%s variables not listed (use --printLinesMax to print more)' % (len(vl)-mx) )
1709
Note: See TracBrowser for help on using the repository browser.