docs: def kippenhahn(self, num_frame, xax, t0_model=0,
                     title='Kippenhahn diagram', tp_agb=0., t_eps=5.e2,
                     plot_star_mass=True, symbol_size=8, c12_bm=False,
                     print_legend=True):
          """Kippenhahn plot as a function of time or model.
  
          Parameters
          ----------
          num_frame : integer
              Number of frame to plot this plot into, if <0 open no new
              figure.
          xax : string
              Either 'model', 'time' or 'logtimerev' to indicate what is
              to be used on the x-axis.
          t0_model : integer, optional
              If xax = 'time' then model for the zero point in time, for
              AGB plots this would be usually the model of the 1st TP,
              which can be found with the Kippenhahn plot.  The default
              is 0.
          title : string, optional
              The figure title.  The default is "Kippenhahn diagram".
          tp_agb : float, optional
              If > 0. then,
              ylim=[h1_min*1.-tp_agb/100 : h1_max*1.+tp_agb/100] with
              h1_min, h1_max the min and max H-free core mass coordinate.
              The default is 0. .
          t_eps : float, optional
              Final time for logtimerev.  The default is '5.e2'.
          plot_star_mass : boolean, optional
              If True, then plot the stellar mass as a line as well.  The
              default is True.
          symbol_size : integer, optional
              Size of convection boundary marker.  The default is 8.
          c12_bm : boolean, optional
              If we plot c12_boundary_mass or not.  The default is False.
          print_legend : boolean, optionla
              Show or do not show legend.  The defalut is True.
  
          """
  
          if num_frame >= 0:
              pyl.figure(num_frame)
  
          t0_mod=[]
  
          if xax == 'time':
              xaxisarray = self.get('star_age')
              if t0_model > 0:
                  ind=self.get('model_number')
                  t0_model=where(ind>t0_model)[0][0]
                  t0_mod=xaxisarray[t0_model]
              else:
                  t0_mod = 0.
              print('zero time is '+str(t0_mod))
          elif xax == 'model':
              xaxisarray = self.get('model_number')
              #t0_mod=xaxisarray[t0_model]
              t0_mod = 0.
          elif xax == 'logtimerev':
              xaxi    = self.get('star_age')
              xaxisarray = np.log10(np.max(xaxi)+t_eps-xaxi)
              t0_mod = 0.
          else:
              print('kippenhahn_error: invalid string for x-axis selction.'+\
                    ' needs to be "time" or "model"')
  
  
          plot_bounds=True
          try:
              h1_boundary_mass  = self.get('h1_boundary_mass')
              he4_boundary_mass = self.get('he4_boundary_mass')
              if c12_bm:
                  c12_boundary_mass = self.get('c12_boundary_mass')
          except:
              try:
                  h1_boundary_mass  = self.get('he_core_mass')
                  he4_boundary_mass = self.get('c_core_mass')
                  if c12_bm:
                      c12_boundary_mass = self.get('o_core_mass')
              except:
                  plot_bounds=False
  
          star_mass         = self.get('star_mass')
          mx1_bot           = self.get('mx1_bot')*star_mass
          mx1_top           = self.get('mx1_top')*star_mass
          mx2_bot           = self.get('mx2_bot')*star_mass
          mx2_top           = self.get('mx2_top')*star_mass
  
  
          if xax == 'time':
              if t0_model>0:
                  pyl.xlabel('$t - t_0$ $\mathrm{[yr]}$')
              else:
                  pyl.xlabel('t / yrs')
          elif xax == 'model':
              pyl.xlabel('model number')
          elif xax == 'logtimerev':
              pyl.xlabel('$\log(t_{final} - t)$  $\mathrm{[yr]}$')
  
          pyl.plot(xaxisarray[t0_model:]-t0_mod,mx1_bot[t0_model:],linestyle='None',color='blue',alpha=0.3,marker='o',markersize=symbol_size,label='convection zones')
          pyl.plot(xaxisarray[t0_model:]-t0_mod,mx1_top[t0_model:],linestyle='None',color='blue',alpha=0.3,marker='o',markersize=symbol_size)
          pyl.plot(xaxisarray[t0_model:]-t0_mod,mx2_bot[t0_model:],linestyle='None',color='blue',alpha=0.3,marker='o',markersize=symbol_size)
          pyl.plot(xaxisarray[t0_model:]-t0_mod,mx2_top[t0_model:],linestyle='None',color='blue',alpha=0.3,marker='o',markersize=symbol_size)
  
          if plot_bounds:
              pyl.plot(xaxisarray[t0_model:]-t0_mod,h1_boundary_mass[t0_model:],color='red',linewidth=2,label='H-free core')
              pyl.plot(xaxisarray[t0_model:]-t0_mod,he4_boundary_mass[t0_model:],color='green',linewidth=2,linestyle='dashed',label='He-free core')
              if c12_bm:
                  pyl.plot(xaxisarray[t0_model:]-t0_mod,c12_boundary_mass[t0_model:],color='purple',linewidth=2,linestyle='dotted',label='C-free core')
          if plot_star_mass is True:
              pyl.plot(xaxisarray[t0_model:]-t0_mod,star_mass[t0_model:],label='$M_\star$')
          pyl.ylabel('$m_\mathrm{r}/\mathrm{M}_\odot$')
          if print_legend:
              pyl.legend(loc=2)
          if tp_agb > 0.:
              h1_min = min(h1_boundary_mass[t0_model:])
              h1_max = max(h1_boundary_mass[t0_model:])
              h1_min = h1_min*(1.-old_div(tp_agb,100.))
              h1_max = h1_max*(1.+old_div(tp_agb,100.))
              print('setting ylim to zoom in on H-burning:',h1_min,h1_max)
              pyl.ylim(h1_min,h1_max)