docs: def _kip(self, cycle_end, mix_thresh, xaxis, sparse):
          """
          *** Should be used with care, therefore has been flagged as
          a private routine ***
          This function uses a threshold diffusion coefficient, above
          which the the shell is considered to be convective, to plot a
          Kippenhahn diagram.
  
          Parameters
          ----------
          cycle_end : integer
              The final cycle number.
          mix_thresh : float
              The threshold diffusion coefficient.
          xaxis : string
              Choose one of 'age', 'cycle', 'log_age' or 'log_time_left'.
          sparse : integer
              Sparsity factor when plotting from cyclelist.
  
          Examples
          --------
          >>> pt=mp.se('/ngpod1/swj/see/mppnp_out/scratch_data/M25.0Z1e-02','.h5')
          >>> pt.kip(10000,'log_time_left',100)
  
          """
  
          original_cyclelist = self.se.cycles
          cyclelist = original_cyclelist[0:cycle_end:sparse]
  
          xx = self.se.ages[:cycle_end:sparse]
          totalmass = []
          m_ini = float(self.se.get('mini'))
  
  
          fig = pl.figure(1)
          ax = pl.subplot(1,1,1)
          fsize = 12
  
          def getlims(d_coeff, massco):
              """
              This function returns the convective boundaries for a cycle,
              given the cycle's dcoeff and massco columns, taking into
              account whether surface or centre are at the top.
  
              """
              plotlims = []
              if massco[0] > massco[-1]:
                  for j in range(-1,-len(d_coeff)-1,-1):
                      if j == -1:
                          if d_coeff[j] >= mix_thresh:
                              plotlims.append(massco[j])
                          else:
                              pass
                      elif (d_coeff[j]-mix_thresh)*(d_coeff[j+1]-mix_thresh) < 0:
                          plotlims.append(massco[j])
                      if j == -len(d_coeff):
                          if d_coeff[j] >= mix_thresh:
                              plotlims.append(massco[j])
                  return plotlims
              else:
                  for j in range(len(d_coeff)):
                      if j == 0:
                          if d_coeff[j] >= mix_thresh:
                              plotlims.append(massco[j])
                          else:
                              pass
                      elif (d_coeff[j]-mix_thresh)*(d_coeff[j-1]-mix_thresh) < 0:
                          plotlims.append(massco[j])
                      if j == len(d_coeff)-1:
                          if d_coeff[j] >= mix_thresh:
                              plotlims.append(massco[j])
                  return plotlims
  
  
          if xaxis == 'age':
              ax.set_xlabel('Age [yr]',fontsize=fsize)
          elif xaxis == 'cycle':
              xx = cyclelist
              ax.set_xlabel('Cycle',fontsize=fsize)
          elif xaxis == 'log_age':
              for i in range(len(xx)):
                  xx[i] = np.log10(xx[i])
              ax.set_xlabel('log$_{10}$(age) [yr]',fontsize=fsize)
          elif xaxis == 'log_time_left':
              for i in range(len(xx)):
                  xx[i] = np.log10(max(xx)-xx[i])
              xx[-2] = xx[-3]-abs(xx[-4]-xx[-3])
              xx[-1] = xx[-2]-abs(xx[-3]-xx[-2])
              ax.set_xlabel('log$_{10}$(time until collapse) [yr]',fontsize=fsize)
  
          #centre-surface flag:
          flag = False
  
          if self.se.get(cyclelist[1],'mass')[0] > self.se.get(cyclelist[1],'mass')[-1]:
              flag = True
  
          for i in range(len(cyclelist)):
              if flag == True:
                  totalmass.append(self.se.get(cyclelist[i],'mass')[0])
              else:
                  totalmass.append(self.se.get(cyclelist[i],'mass')[-1])
              percent = int(i*100/len(cyclelist))
              sys.stdout.flush()
              sys.stdout.write("\rcreating color map " + "...%d%%" % percent)
              d_coeff = self.se.get(cyclelist[i],'dcoeff')
              massco = self.se.get(cyclelist[i],'mass')
              plotlims = getlims(d_coeff,massco)
              for k in range(0,len(plotlims),2):
                  ax.axvline(xx[i],ymin=old_div(plotlims[k],m_ini),ymax=old_div(plotlims[k+1],m_ini),color='b',linewidth=0.5)
  
  
          ax.plot(xx, totalmass, color='black', linewidth=1)
          if xaxis == 'log_time_left':
              ax.axis([xx[0],xx[-1],0.,m_ini])
          else:
              ax.axis([min(xx),max(xx),0.,m_ini])
          ax.set_ylabel('Mass [$M_{\odot}$]',fontsize=fsize)
  
          pl.show()