docs: def plot_isoratios(self,xiso,yiso,fign=1,spec=None,deltax=True,deltay=True,logx=False,logy=False,
                          title=None,legend=None,legloc='lower right',errbar=True,dcycle=500,addiso=None,
                          co_toggle='c',cust_toggle=None,shift=0,weighting=None,zoneselect=None,iniabufile='iniab2.0E-02GN93.ppn',
                          plt_sparse=1,plt_symb='o',plt_col='k',plt_lt='-',plt_lw=1.,alpha_dum=1.,plt_massrange=False,plt_show=True,
                          figsave=False):
          '''
          This is the new routine to plot isotopic ratios for ALL input. rt, June 2014
  
          Parameters:
          -----------
          xiso : np.array
              x data to plot. This can be an array or a list of arrays, depending on who calls the routine
          yiso : np.array
              y data to plot. This can be an array or a list of arrays, depending on who calls the routine
          fign : integer, optional
              Figure number
          spec : string, optional
              What specifications do you want to do when coming from nugridse models. Choose 'surf' for
              surface models or 'exp' for explosions (out files)
          deltax : boolean, optional
              X axis in delta values?
          deltay : boolean, optional
              Y axis in delta values?
          logx : boolean, optional
              Logarithmic x axis?
          logy : boolean, optional
              Logarithmic y axis?
          title : string, optional
              Title for your plot
          legend : string, optional
              Legend for your model / grains. For grains the legend is automatically taken from the
              grain class
          legloc : string / integer, optional
              Location of the legend, use matplotlib standard. Use None to not plot legend if plotted
              by default, e.g., from grain class routine.
          errbar : boolean, optional
              Error bars on grain data?
          dcycle : integer, optional
              Difference between cycles to take for thermal pulse searching, if searching is
              deactivated, dcycle describes how often cycles are sampled.  The default is 500.
          addiso : list, optional
              For explosive models. Add an isotope.  Format ['C-12', 0.5 ,'N-12'] to add N12
              to C12 and multiply it with a factor of 0.5.  Multiple isotopes can be added, the
              factor is optional and does not have to be given.  Isotopes can be added to other
              isotopes as well, i.e., [['C-12', 'N-12'], ['C-13', 'N-13']].  The default None.
              Notice that while addiso = [['N-14','O-14'],['N-14',fractionation,'C-14']]
              works, other options like addiso = [['N-14','O-14',fractionation,'C-14']] or
              addiso = [['N-14',1,'O-14',fractionation,'C-14']] are not working and give Typerror.
              CAREFUL, that for the option addiso = [['N-14',fractionation,'C-14','O-14']] there is
              no error message, but the fractionation is applied to both O14 and C14!
          co_toggle : string, optional
              For explosive models, choose what shells you want to look for! Select 'c' for
              selecting zones with C/O >= 1.  Select 'o' for C/O <= 1.  If 'a' takes the
              whole star.  The defalut is 'c'. See cust_toggle (below) for an alternative!
          cust_toggle : list, optional
              This option is like co_toggle (and overwrites it when chosen) but lets you choose
              your own comparison. For example you want to find zones that have a 10 fold
              overabundance of Ti-46 and Ti-47 over O-16 and Zr-96, you can choose here
              [['Ti-46','Ti-47'],['O-16','Zr-96'],100.] Assuming the first list is is x, the
              second list y, and the comparator number is f, the statement only plots shells
              in which the condition x/y>f is fulfilled. x and y are number sums of the chosen
              isotopes, f has to be given as a float. This is only for explosive shells. Please
              note, if this toggle is NOT None, then co_toggle is overwritten!
          shift : integer, optional
              For explosive models, how much do you want to shift the models back from the
              last cycle? By default (0) the last cycle is taken.
          weighting : string, optional
              For explosive models. If None then, plot every profile separately.  If 'zone'
              then, average each zone.  If 'all' then average all selected zones.  The
              default is None.
          zoneselect: string, optional
              For explosive models. Select if you want to plot 'all' zones or outer most zone.
              Arguments are 'all' and 'top', respectively. Default is None, then the user is
              asked to provide this information during the routine as input.
          iniabufile : string, optional
              Initial abundance file. Use absolute path for your file or filename to choose a
              given file in USEEPP. Attention: You need a standard tree checked out from SVN
          plt_sparse : integer, optional
              Every how many datapoints is the plot done? Not used for some routines!
          plt_symb : string, optional
              Symbol for the plot. In case of grains, this is handled automatically.
          plt_col : string / float, optional
              Color for plotted curve. In case of grains, this is handled automatically.
          plt_lt : string, optional
              line type for plot.
          plt_lw : float, optional
              Line width for plot.
          alpha_dum : trasparency to apply to grains data, in case of many data are plotted.
              This may be allpied also for theoretical curves.
          plt_massrange : boolean, optional
              For explosive models. Plot mass of shell with first and last datapoint of
              each zone.  If list given, label those zones.  The default is False.
          plt_show : boolean, optional
              Do you want to show the plot or not?
          figsave : string, optional
              Give path and filename here, if you want to save the figure.
          '''
  
          from . import utils as u
  
          ### WORK ON PATH ###
          # define svn path form path where script runs, depending on standard input or not
          if len(iniabufile.split('/')) == 1 :   # means not an absolute path!
              scriptpathtmp = __file__
              if len(scriptpathtmp.split('/')) == 1:   # in folder where nugridse is
                  scriptpathtmp = os.path.abspath('.') + '/nugridse.py'   # to get the current dir
              svnpathtmp = '/'
              for i in range(len(scriptpathtmp.split('/'))-3):   # -3 to go to folders up!
                  if scriptpathtmp.split('/')[i] != '':
                      svnpathtmp += scriptpathtmp.split('/')[i] + '/'
              iniabufile = svnpathtmp + 'frames/mppnp/USEEPP/' + iniabufile   # make absolute path for iniabufile
  
  
          ### get solar system ratios for the isotopes that are specified in the input file ###
          inut = u.iniabu(iniabufile)
          try:
              xrat_solsys = inut.isoratio_init(xiso)
          except KeyError:   # if isotope not available, e.g., if plotting Ti-44 / Ti-48 ratio
              xrat_solsys = 0.
          try:
              yrat_solsys = inut.isoratio_init(yiso)
          except KeyError:
              yrat_solsys = 0.
  
          # number ratio for solar system ratio
          xrat_solsys *= (old_div(float(xiso[1].split('-')[1]), float(xiso[0].split('-')[1])))
          yrat_solsys *= (old_div(float(yiso[1].split('-')[1]), float(yiso[0].split('-')[1])))
  
  
          # initialize xdataerr and ydataerr as None
          xdataerr = None
          ydataerr = None
  
  
  
          ### DO PLOTS FROM NUGRIDSE CLASS ###
          if self._classTest() == 'se':
              if spec==None:
                  spec = str(eval(input('Please specify \'surf\' for surface models (AGB stars) or \'exp\' for explosive'
                               'models and zone finding, etc., and press enter: ')))
  
              ### SURFACE MODELS - PLOT AGB STAR STUFF ###
              if spec == 'surf':
                  print('Plotting AGB star stuff')
  
                  # read in thermal pulse position and co ratio
                  tp_pos, co_return = self._tp_finder(dcycle)
                  tp_pos_tmp = []
                  co_return_tmp  = []
                  tp_pos_tmp.append(1)
                  co_return_tmp.append(co_return[0])
                  for i in range(len(tp_pos)):
                      tp_pos_tmp.append(tp_pos[i])
                      co_return_tmp.append(co_return[i])
                  tp_pos = tp_pos_tmp
                  co_return = co_return_tmp
  
                  # read in data
                  iso_alldata = self.get(tp_pos,[xiso[0],xiso[1],yiso[0],yiso[1]])
                  xrat = np.zeros(len(iso_alldata))
                  yrat = np.zeros(len(iso_alldata))
                  for i in range(len(iso_alldata)):
                      xrat[i] = old_div(iso_alldata[i][0], iso_alldata[i][1])
                      yrat[i] = old_div(iso_alldata[i][2], iso_alldata[i][3])
  
                  # make number ratios
                  for i in range(len(xrat)):
                      xrat[i] *= old_div(float(xiso[1].split('-')[1]), float(xiso[0].split('-')[1]))
                      yrat[i] *= old_div(float(yiso[1].split('-')[1]), float(yiso[0].split('-')[1]))
  
                  # if delta values are requested, need to calculate those now
                  if deltax:
                      xrat = (old_div(xrat,xrat_solsys) -1.) * 1000.
                  if deltay:
                      yrat = (old_div(yrat,yrat_solsys) -1.) * 1000.
  
                  # now we might have o and c rich zones. prepare stuff for plotting
                  xdata_o = []
                  ydata_o = []
                  xdata_c = []
                  ydata_c = []
                  for i in range(len(co_return)):
                      if co_return[i] <= 1.:
                          xdata_o.append(xrat[i])
                          ydata_o.append(yrat[i])
                      else:
                          xdata_c.append(xrat[i])
                          ydata_c.append(yrat[i])
                  if xdata_o != [] and xdata_c != []:
                      xdata_o.append(xdata_c[0])
                      ydata_o.append(ydata_c[0])
                  # now make the styles
                  style_o = [plt_symb + '--', plt_col, '1.', '4', '2',None]
                  style_c = [plt_symb + plt_lt,  plt_col, plt_col, '7.', '1', legend]
  
                  # now make data for plotting
                  xdata = []
                  ydata = []
                  style = []
                  if xdata_o != []:
                      xdata.append(xdata_o)
                      ydata.append(ydata_o)
                      style.append(style_o)
                  if xdata_c != []:
                      xdata.append(xdata_c)
                      ydata.append(ydata_c)
                      style.append(style_c)
  
  
  
  
              ### EXPLOSIVE MODELS ###
              elif spec == 'exp':
                  print('explosive models')
                  # compatibility
                  co_toggle = co_toggle.lower()
                  isotope_list = [xiso[0],xiso[1],yiso[0],yiso[1]]
  
                  # cycle
                  cyc_no = self.se.cycles[len(self.se.cycles)-1-shift]
                  mco_data = self.get(cyc_no,['mass','C-12','C-13','O-16','O-17','O-18',xiso[0],xiso[1],yiso[0],yiso[1]])
                  mass = mco_data[0]
                  # if no custom toggle for enrichment
                  if cust_toggle == None:
                      c_elem = mco_data[1]+mco_data[2]
                      o_elem = mco_data[3]+mco_data[4]+mco_data[5]
                      co_ratio = c_elem / o_elem * (old_div(16., 12.))
                      co_comp_val = 1.
                  else:
                      co_data1 = self.get(cyc_no,cust_toggle[0])
                      co_data2 = self.get(cyc_no,cust_toggle[1])
                      for i in range(len(co_data1)):
                          if i == 0:
                              c_elem = co_data1[i]
                          else:
                              c_elem += co_data1[i]
                      for i in range(len(co_data2)):
                          if i == 0:
                              o_elem = co_data2[i]
                          else:
                              o_elem += co_data2[i]
                      # now we need to make the mass number of everything in here to make number ratios
                      massn1 = 0.
                      for i in range(len(co_data1)):
                          massn1 += sum(co_data1[i]) * float(cust_toggle[0][i].split('-')[1])
                      massn1 /= sum(c_elem)
                      massn2 = 0.
                      for i in range(len(co_data2)):
                          massn2 += sum(co_data2[i]) * float(cust_toggle[1][i].split('-')[1])
                      massn2 /= sum(o_elem)
                      co_ratio = c_elem / o_elem * (old_div(massn2, massn1))   # this has nothing to do with a C/O ratio anymore! but keep name
                      # comparator value
                      co_comp_val = float(cust_toggle[2])
  
                  # get the data now
                  isotope_profile = []
                  for i in range(6,10):   # in mco_data
                      isotope_profile.append(mco_data[i])
  
                  # add radioactive isotopes (if given)
                  if addiso != None:
                      if type(addiso[0] == list):   # then list of lists
                          for i in range(len(addiso)):
                              for j in range(4):
                                  if isotope_list[j] == addiso[i][0]:
                                      multiplicator_addiso = 1.
                                      try:
                                          multiplicator_addiso = float(addiso[i][1])
                                          starter = 2
                                      except ValueError:
                                          starter = 1
                                      for k in range(starter,len(addiso[i])):
                                          isotope_profile[j] += array(self.get(cyc_no,addiso[i][k])) * multiplicator_addiso
                      else:
                          for j in range(4):
                              if isotope_list[j] == addiso[0]:
                                  multiplicator_addiso = 1.
                                  try:
                                      multiplicator_addiso = float(addiso[1])
                                      starter = 2
                                  except ValueError:
                                      starter = 1
                                  for k in range(starter,len(addiso)):
                                      isotope_profile[j] += array(self.get(cyc_no,addiso[k])) * multiplicator_addiso
  
                  # search for carbon / oxygen rich layers
                  crich = []   # alternating start stop values. if odd number, then surface is c-rich, but add stop number
                  dumb = True
                  if cust_toggle != None:
                      for i in range(len(co_ratio)):
                          if dumb:
                              if co_ratio[i] >= co_comp_val:
                                  crich.append(i)
                                  dumb = False
                                  continue
                          else:
                              if co_ratio[i] < co_comp_val:
                                  crich.append(i)
                                  dumb = True
                  elif co_toggle != 'a':
                      for i in range(len(co_ratio)):
                          if co_toggle == 'c':   # carbon rich
                              if dumb:
                                  if co_ratio[i] >= 1:
                                      crich.append(i)
                                      dumb = False
                                      continue
                              else:
                                  if co_ratio[i] < 1:
                                      crich.append(i)
                                      dumb = True
                          elif co_toggle == 'o':   # oxygen rich
                              if dumb:
                                  if co_ratio[i] <= 1:
                                      crich.append(i)
                                      dumb = False
                                      continue
                              else:
                                  if co_ratio[i] > 1:
                                      crich.append(i)
                                      dumb = True
                          else:
                              print('Select your enrichment!')
                              return None
                  else:   # take whole star
                      print('Using all profiles to mix')
                      crich.append(0)
                      crich.append(len(co_toggle))
  
                  if len(crich)%2 == 1:
                      crich.append(len(co_ratio)-1)
  
                  if len(crich) == 0:
                      print('Star did not get rich in C or O, depending on what you specified')
                      return None
  
                  # make isotope_profile into array and transpose
                  isotope_profile = array(isotope_profile).transpose()
  
  
                  # Ask user which zones to use
                  if co_toggle != 'a':
                      if cust_toggle != None:
                          print('\n\nI have found the following zones:\n')
                      elif co_toggle == 'c':
                          print('\n\nI have found the following carbon rich zones:\n')
                      elif co_toggle == 'o':
                          print('\n\nI have found the following oxygen rich zones:\n')
  
                      mass_tmp = zeros((len(crich)))
                      for i in range(len(crich)):
                          mass_tmp[i] = mass[crich[i]]
  
                      j = 1
                      for i in range(old_div(len(crich),2)):
                          print('Mass range (' + str(j) + '):\t' + str(mass_tmp[2*i]) + ' - ' + str(mass_tmp[2*i+1]))
                          j += 1
                      print('\n')
                      if zoneselect == 'all':
                          usr_zones = 0
                      elif zoneselect == 'top':
                          usr_zones = [j-1]
                      else:
                          usr_zones = eval(input('Please select which mass range you want to use. Select 0 for all zones. Otherwise give one zone or a list of zones separated by comma (e.g.: 1, 2, 4): '))
  
                      crich_dumb = crich
                      if usr_zones == 0:
                          print('I continue w/ all zones then')
                      elif type(usr_zones) == int:   # only one zone selected
                          tmp = int(usr_zones)-1
                          crich = crich_dumb[2*tmp:2*tmp+2]
                      else:
                          crich = []
                          for i in range(len(usr_zones)):
                              tmp = int(usr_zones[i])-1
                              crich.append(crich_dumb[2*tmp])
                              crich.append(crich_dumb[2*tmp + 1])
  
                  # weight profiles according to weighting factor using the selected crich
                  # define isos_to_use variable for later
                  if weighting == None:
                      isos_to_use = []
                      for i in range(old_div(len(crich),2)):
                          isos_dumb = []
                          n = crich[2*i]
                          while n <= crich[2*i+1]:
                              isos_dumb.append(isotope_profile[n])
                              n += 1
                          isos_to_use.append(array(isos_dumb))
  
  
                  elif weighting.lower() == 'zone' or weighting.lower() == 'zones':
                      # make array w/ mass weigted isotope ratio (4) for all mass zones
                      isotope_profile_cweight = zeros((old_div(len(crich),2),4))
                      mass_tot = []
                      for i in range(len(isotope_profile_cweight)):   # 2*i is start, 2*i+1 is stop value
                          if crich[2*i] == 0:
                              print('C- / O-rich in first shell (core).')
                          else:
                              dumb = crich[2*i + 1]
                              j = crich[2*i]
                              mass_tmp = 0
                              while j <= dumb:
                                  mass_shell = mass[j] - mass[j-1]
                                  mass_tmp += mass_shell
                                  for k in range(4):
                                      isotope_profile_cweight[i][k] += isotope_profile[j][k]*mass_shell
                                  j += 1
                              mass_tot.append(mass_tmp)
                      for i in range(len(isotope_profile_cweight)):
                          for j in range(4):
                              isotope_profile_cweight[i][j] /= mass_tot[i]
                      isos_to_use = [array(isotope_profile_cweight)]
  
                  elif weighting.lower() == 'all':   # average all zones by mass
                      isos_tmp = zeros((1, len(isotope_profile[0])))
                      for i in range(len(isotope_profile)-1):   # neglect surface effects
                          for j in range(len(isos_tmp[0])):
                              mass_shell = mass[i+1] - mass[i]
                              isos_tmp[0][j] += isotope_profile[i][j]*mass_shell
                      # weight all
                      isos_tmp /= sum(mass)
                      isos_to_use = [isos_tmp]
  
  
                  # change to isotope numbers from mass!
                  for i in range(len(isos_to_use)):
                      for j in range(len(isos_to_use[i])):
                          for k in range(len(isos_to_use[i][j])):
                              # here we just divide 'iso_massf' output with the mass number
                              # this means that in the end, the isotope ratios in number space are correc
                              # but have to use ratios from here on for meaningful stuff
                              isos_to_use[i][j][k] /= float(isotope_list[k].split('-')[1])
  
                  # do the ratios and stuff
                  ratiox = []
                  ratioy = []
                  for i in range(len(isos_to_use)):
                      ratiox_dumb = []
                      ratioy_dumb = []
                      for j in range(len(isos_to_use[i])):
                          ratiox_dumb.append(old_div(isos_to_use[i][j][0], isos_to_use[i][j][1]))
                          ratioy_dumb.append(old_div(isos_to_use[i][j][2], isos_to_use[i][j][3]))
                      ratiox.append(array(ratiox_dumb))
                      ratioy.append(array(ratioy_dumb))
  
                  # make arrays for ratiox and ratioy
                  ratiox = array(ratiox)
                  ratioy = array(ratioy)
  
                  # make number ratio out of everything
                  for i in range(len(ratiox)):
                      for j in range(len(ratiox[i])):
                          ratiox[i][j] *= (old_div(float(xiso[1].split('-')[1]), float(xiso[0].split('-')[1])))
                          ratioy[i][j] *= (old_div(float(yiso[1].split('-')[1]), float(yiso[0].split('-')[1])))
  
                  if deltax:
                      ratiox_tmp = []
                      for i in range(len(ratiox)):
                          ratiox_tmp_tmp = []
                          for j in range(len(ratiox[i])):
                              ratiox_tmp_tmp.append((old_div(ratiox[i][j], xrat_solsys) - 1.) * 1000.)
                          ratiox_tmp.append(ratiox_tmp_tmp)
                      ratiox = array(ratiox_tmp)
                  if deltay:
                      ratioy_tmp = []
                      for i in range(len(ratioy)):
                          ratioy_tmp_tmp = []
                          for j in range(len(ratioy[i])):
                              ratioy_tmp_tmp.append((old_div(ratioy[i][j], yrat_solsys) - 1.) * 1000.)
                          ratioy_tmp.append(ratioy_tmp_tmp)
                      ratioy = array(ratioy_tmp)
  
                  # create massrange array if necessary
                  plt_massrange_lst = []
                  if plt_massrange==True:   # use == True because otherwise the list enters here too... why?
                      for i in range(len(ratiox)):
                          plt_massrange_lst.append([ratiox[i][0], ratioy[i][0], mass[crich[2*i]]])
                          plt_massrange_lst.append([ratiox[i][len(ratiox[i])-1], ratioy[i][len(ratioy[i])-1], mass[crich[2*i+1]]])   # start: x-ratio, y-ratio, mass label
  
                  elif plt_massrange != False:
                      for plt_mr_val in plt_massrange:
                          mrng_i = 0
                          while plt_mr_val > mass[mrng_i] and mrng_i < len(mass):
                              mrng_i += 1
                          if mrng_i > 0:
                              mrng_i -= 1
                          mratx = old_div(isotope_profile[mrng_i][0],isotope_profile[mrng_i][1])
                          mraty = old_div(isotope_profile[mrng_i][2],isotope_profile[mrng_i][3])
                          if deltax:
                              mratx = (old_div(mratx, ratiox_solsys) - 1.) * 1000.
                          if deltay:
                              mraty = (old_div(mraty, ratioy_solsys) - 1.) * 1000.
                          plt_massrange_lst.append([mratx,mraty,mass[mrng_i]])
  
                  # make style and prepare for plotting here
                  xdata = ratiox
                  ydata = ratioy
                  style_tmp0= [plt_symb + plt_lt,  plt_col, plt_col, '13.', '1', legend]
                  style_tmp = [plt_symb + plt_lt,  plt_col, plt_col, '13.', '1', None]
                  style = []
                  for i in range(len(xdata)):
                      if i == 0:
                          style.append(style_tmp0)
                      else:
                          style.append(style_tmp)
  
              else:
                  print('You did not specify a useful spec argument -> abort.')
                  return None
  
          ### PLOTS FROM GRAIN CLASS ###
          if self._classTest() == 'grain':
              print('Presolar grains are cool!')
              xdata,xdataerr,ydata,ydataerr,style = self.plot_ratio_return(xiso,yiso,deltax,deltay)
              legend=True
              plt_sparse=1   # to avoid monkey input
              plt_lw = 0.
  
          ### PLOT ###
          # data is prepared now, make the plots. data must be in format
          # [[data1],[data2],[data3],...]
          # three arrays like this, for xdata, ydata, and style.
          # style format: symbol, edge color, face color, symbol size, edge width, label
          # this is then compatible with grain.py style definitions
  
          # Size of font etc.
          params = {'axes.labelsize':  20,
                    'text.fontsize':   14,
                    'legend.fontsize': 14,
                    'xtick.labelsize': 14,
                    'ytick.labelsize': 14}
          pl.rcParams.update(params)
  
          pl.figure(fign)
          for i in range(len(xdata)):
              if errbar:
                  if xdataerr != None or ydataerr != None:
                      pl.errorbar(xdata[i],ydata[i],xerr=xdataerr[i],yerr=ydataerr[i],marker=style[i][0],
                                  color=style[i][1],linestyle='',lw=2,markevery=plt_sparse,alpha=alpha_dum)
              pl.plot(xdata[i],ydata[i],style[i][0],c=style[i][1],mfc=style[i][2],ms=float(style[i][3]),
                      mew=float(style[i][4]),label=style[i][5],markevery=plt_sparse,linewidth=plt_lw,alpha=alpha_dum)
  
          # plot text labels if necessary
          if plt_massrange != False:
              for mrng_ind in range(len(plt_massrange_lst)):
                  pl.text(plt_massrange_lst[mrng_ind][0], plt_massrange_lst[mrng_ind][1],
                          str(round(plt_massrange_lst[mrng_ind][2], 2)), ha='right', va='bottom', color=plt_col,fontsize=15.)
  
  
          # log?
          if logx and logy == False:
              pl.semilogx()
          elif logx == False and logy:
              pl.semilogy()
          elif logx and logy:
              pl.loglog()
  
          # legend
          if legend != None and legloc != None:
              pl.legend(loc=legloc)
  
          # title and labels
          if title != None:
              pl.title(title)
          if deltax:
              pl.xlabel('$\delta$($^{' + xiso[0].split('-')[1] + '}$' +xiso[0].split('-')[0] + '/$^{' + xiso[1].split('-')[1] + '}$' +xiso[1].split('-')[0] + ')' )
          else:
              pl.xlabel('$^{' + xiso[0].split('-')[1] + '}$' +xiso[0].split('-')[0] + '/$^{' + xiso[1].split('-')[1] + '}$' +xiso[1].split('-')[0])
          if deltay:
              pl.ylabel('$\delta$($^{' + yiso[0].split('-')[1] + '}$' +yiso[0].split('-')[0] + '/$^{' + yiso[1].split('-')[1] + '}$' +yiso[1].split('-')[0] + ')' )
          else:
              pl.ylabel('$^{' + yiso[0].split('-')[1] + '}$' +yiso[0].split('-')[0] + '/$^{' + yiso[1].split('-')[1] + '}$' +yiso[1].split('-')[0])
  
          # plot horizontal and vertical lines
          print(xrat_solsys, yrat_solsys)
          if deltay:
              pl.axhline(0,color='k')
          else:
              pl.axhline(yrat_solsys,color='k')
          if deltax:
              pl.axvline(0,color='k')
          else:
              pl.axvline(xrat_solsys,color='k')
  
          # borders of plot
          pl.gcf().subplots_adjust(bottom=0.15)
          pl.gcf().subplots_adjust(left=0.15)
  
          # save and show
          if figsave != False:
              pl.savefig(figsave)
  
          if plt_show:
              pl.show()