#!/usr/bin/env python
# coding: utf-8

# *Examples are given here in abbreviated form. We invite you to look over our other examples which explain these calculations in greater detail, as well as read the freud documentation to better understand these examples*

# In[1]:


from aps_helpers import *
from freud import parallel
parallel.setNumThreads(4)


# # Visualization
# 
# While *Freud* does not provide any methods for visualization, *Freud* allows you to enhance your visualizations in ways which make them easier to understand. Here we demonstrate how to visualize a trajectory for a simple demonstration system.
# 
# ### System Information
# 
# - Simulated with *HOOMD-Blue*
#     - Hard Particle Monte Carlo (hpmc) regime
# - \(N = 4096\) particles
# - *NVT* Thermodynamic ensemble
# - Example data contains 400 trajectory frames
# 
# ## Load the data
# 
# *Freud* does not contain any file readers or writers. For this demonstration I have saved the relevant information as *numpy binary* files.

# ## Visualize the data
# 
# Now we visualize our system, using bokeh to draw our shapes

# In[2]:


# load the data
# data_path = "ex_data/phi065"
# fluid_data = DemoData(data_path)
# fluid_viz = fluid_data.plot_frame(-1, title="65% Density")

data_path = "ex_data/phi075"
fluid_data = DemoData(data_path)
fluid_viz = fluid_data.plot_frame(-1, title="75% Density")
# create a grid for this

data_path = "ex_data/phi085"
solid_data = DemoData(data_path)
solid_viz = solid_data.plot_frame(-1, title="85% Density", linked_plot=fluid_viz)

r = row(fluid_viz, solid_viz)
show(r)


# ## Nearest Neighbors and Hexatic Calculation
# 
# Just viewing particle positions and orientations doesn't give us much information. Let's consider two different ways of looking at our system in more detail: the number of nearest neighbors a particle has, and the angle of the hexatic order parameter for a particle.
# 
# ### Nearest Neighbors Calculation
# 
# First, let's view the number of nearest neighbors each particle has.

# In[3]:


# import freud locality object
from freud import locality
# set number of neighbors
n_neigh = 6
# create freud nearest neighbors object
nn = locality.NearestNeighbors(rmax=1.5,
    n_neigh=n_neigh,strict_cut=True)

# compute nearest neighbors for 6 nearest neighbors
nn.compute(fluid_data.freud_box(-1),
           fluid_data.pos_data[-1],
           fluid_data.pos_data[-1])
# get the computed neighborlist
n_list = np.copy(nn.getNeighborList())
# get the number of particles
num_particles = nn.getNRef()
# display nearest neighbors for particle 1000
neigh_fluid = fluid_data.plot_single_neighbor(-1, 1000, n_list, num_particles, title="75% Density")
# compute nearest neighbors for 6 nearest neighbors
nn.compute(solid_data.freud_box(-1),
           solid_data.pos_data[-1],
           solid_data.pos_data[-1])
# get the computed neighborlist
n_list = np.copy(nn.getNeighborList())
# get the number of particles
num_particles = nn.getNRef()
# display nearest neighbors for particle 1000
neigh_solid = solid_data.plot_single_neighbor(-1, 1000, n_list, num_particles, title="85% Density")
r = row(neigh_fluid, neigh_solid)
show(r)


# In[4]:


# set number of neighbors
n_neigh = 6
# create freud nearest neighbors object
# set rmax to some value
rmax=1.45
nn = locality.NearestNeighbors(rmax=rmax,
    n_neigh=n_neigh,strict_cut=True)

# compute nearest neighbors for 6 nearest neighbors
nn.compute(fluid_data.freud_box(-1),
           fluid_data.pos_data[-1],
           fluid_data.pos_data[-1])
# get the computed neighborlist
n_list = np.copy(nn.getNeighborList())
# get the number of particles
num_particles = nn.getNRef()
# display nearest neighbors for particle 1000
f_neighbor = fluid_data.plot_neighbors(-1, n_list, num_particles, n_neigh, title="75% Density")

# compute nearest neighbors for 6 nearest neighbors
nn.compute(solid_data.freud_box(-1),
           solid_data.pos_data[-1],
           solid_data.pos_data[-1])
# get the computed neighborlist
n_list = np.copy(nn.getNeighborList())
# get the number of particles
num_particles = nn.getNRef()
# display nearest neighbors for particle 1000
s_neighbor = solid_data.plot_neighbors(-1, n_list, num_particles, n_neigh, linked_plot=f_neighbor, title="85% Density")

r = row(f_neighbor, s_neighbor)
show(r)


# ### Hexatic Order Parameter
# 
# In addition to viewing the number of neighbors to find defects, *freud* can also compute a variety of order parameters. In this example we color each hexagon by its hexatic order parameter angle so as to better see the overall orientation of neighbors in the system.

# In[5]:


# import freud order module
from freud import order
# create hexatic order parameter object
hex_order = order.HexOrderParameter(rmax=1.5, k=6, n=6);
# compute hexatic order for 6 nearest neighbors
hex_order.compute(fluid_data.freud_box(-1),
    fluid_data.pos_data[-1])
# get values from freud object
psi_k = hex_order.getPsi()
avg_psi_k = np.mean(psi_k)

h_fluid = fluid_data.plot_hexatic(-1, psi_k, avg_psi_k, linked_plot=f_neighbor, title="75% Density")

# compute hexatic order for 6 nearest neighbors
hex_order.compute(solid_data.freud_box(-1),
    solid_data.pos_data[-1])
# get values from freud object
psi_k = hex_order.getPsi()
avg_psi_k = np.mean(psi_k)

h_solid = solid_data.plot_hexatic(-1, psi_k, avg_psi_k, linked_plot=f_neighbor, title="85% Density")
r = row(h_fluid, h_solid)
show(r)


# ## Compare with visualizing by particle orientation

# In[6]:


a_fluid = fluid_data.plot_orientation(-1, linked_plot=f_neighbor, title="75% Density")
a_solid = solid_data.plot_orientation(-1, linked_plot=f_neighbor, title="85% Density")
r = row(a_fluid, a_solid)
show(r)


# In[7]:


g = gridplot([[h_fluid, a_fluid, f_neighbor],[h_solid, a_solid, s_neighbor]])
show(g)