histo_pca.py 1.1 KB

12345678910111213141516171819202122232425262728293031323334
  1. """Histogram along a PCA axis"""
  2. import numpy as np
  3. from vedo import Points, pca_ellipse, Arrow2D, Goniometer
  4. from vedo.pyplot import Figure, histogram
  5. data = np.random.randn(1000, 3)
  6. pts = Points(data).color('#1f77b4').ps(6)
  7. pts.scale([2,1,0.01]).rotate_z(45).shift(5,1) # rotate and shift!
  8. # Recover the rotation pretending we only know the points
  9. # Fit a 1-sigma ellipse to the points
  10. elli = pca_ellipse(pts)
  11. ec, e1, e2 = elli.center, elli.axis1, elli.axis2
  12. arrow1 = Arrow2D(ec, ec - 3*e1)
  13. arrow2 = Arrow2D(ec, ec + 3*e2)
  14. angle = np.arctan2(e1[1], e1[0]) * 180/np.pi
  15. mypts = pts.clone() # rotate back to make the histo:
  16. mypts.shift(-ec).rotate_z(-angle)
  17. histo = histogram( # a Histogram1D(Figure) object
  18. mypts.points[:,1], # grab the y-values (PCA2)
  19. ytitle='', title=' ', # no automatic title, no y-axis
  20. c='#1f77b4', # color
  21. aspect=16/9, # aspect ratio
  22. )
  23. histo.rotate_z(90 + angle).pos(ec - 6*e1)
  24. gon = Goniometer(ec-5.5*e1, ec, [ec[0]-5.5*e1[0], ec[1],0]).z(0.2)
  25. fig = Figure([0,14], [-4,9], aspect="equal", title=__doc__)
  26. fig += [pts, elli.z(-0.1), arrow1, arrow2, gon, histo]
  27. fig.show(zoom='tight').close()