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

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

01.00.00

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