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

# ChEn-3170: Computational Methods in Chemical Engineering Fall 2018 UMass Lowell; Prof. V. F. de Almeida **16Sep2018**
# 
# # 04. Arrays Operations

# ---
# ## Table of Contents
# * [Vectors](#vectors)
# * [Matrices](#matrices)
# ---

# In[1]:


'''Python packages are accessed with an import directive as such:'''

import numpy as np  # import the package and create the alias: np


# ## Vectors<a id="vectors"></a>

# In[2]:


'''Element-by-element addition or subtraction'''

vec1 = np.array( np.random.random(5) )
print('vec1     =',vec1)
vec2 = np.array(np.random.random(5))
print('vec2     =',vec2)
result = vec1 + vec2         # element-by-element sum
print('addition   =',result)
result = vec1 - vec2         # element-by-element subtraction
print('difference =',result)


# In[3]:


'''Element-by-element product or division'''

vec1 = np.array(np.random.random(5))
print('vec1    =',vec1)
vec2 = np.array(np.random.random(5))
print('vec2    =',vec2)
result = vec1 * vec2        # element-by-element product
print('product  =',result)
result = vec1 / vec2        # element-by-element division
print('division =',result)


# In[4]:


'''Vector scalar product'''

vec1 = np.array(np.random.random(5))
print('vec1           =',vec1)
vec2 = np.array(np.random.random(5))
print('vec2           =',vec2)
result = np.dot(vec1, vec2)     # scalar or dot product
print('scalar product =',result)


# In[7]:


'''Scaling of a vector'''

vec = np.array(np.random.random(5))
print('vec    =',vec)
factor = 0.345
scaled = factor * vec     # scaling of vec element-by-element product
print('scaled =', scaled) # assigned to new variable `scaled`

vec *= factor          # in-place scaling
print('vec    =',vec)


# In[8]:


'''Mathematical Operations on a Vector'''

vec = np.array(np.random.random(5))
print('vec      =',vec)

log_vec = np.log(vec)         # natural log element-by-element
print('log(vec) =',log_vec)

exp_vec = np.exp(log_vec)     # exponential
print('exp(vec) =',exp_vec)

sin_vec = np.sin(vec)         # sine
print('sin(vec) =',sin_vec)

vec_cubed = vec**3            # powers
print('vec^3    =',vec_cubed)

vec_mean = vec.mean()         # arithmetic mean
print('mean(vec) =',vec_mean)

vec_std = vec.std()           # standard deviation
print('std(vec) =',vec_std)


# ## Matrices<a id="matrices"></a>

# In[9]:


'''Element-by-element addition or subtraction'''

mat1 = np.random.random( (3,3) )
print('mat1     =\n',mat1)
mat2 = np.random.random( (3,3) )
print('mat2     =\n',mat2)
result = mat1 + mat2              # element-by-element sum
print('addition   =\n',result)
result = mat1 - mat2              # element-by-element subtraction
print('difference =\n',result)


# In[10]:


'''Element-by-element product or division'''

mat1 = np.random.random((3,3))
print('mat1     =\n',mat1)
mat2 = np.random.random((3,3))
print('mat2     =\n',mat2)
result = mat1 * mat2          # element-by-element product
print('product   =\n',result)
result = mat1 / mat2          # element-by-element division (cross your fingers)
print('division =\n',result)


# In[14]:


'''Produce Noise on a Matrix Image'''

from matplotlib import pyplot as plt     # import the pyplot function of the matplotlib package
plt.rcParams['figure.figsize'] = [20, 4] # extend the figure size on screen output
# copy URL of image in the images/ course repository;
# click on the image link, use the Download buttom to create this URL
block = plt.imread('https://raw.githubusercontent.com/dpploy/chen-3170/master/notebooks/images/glacier.png',format='png')

plt.figure(1)
plt.imshow( block )
plt.title('Matrix Reloaded',fontsize=14)
print('mtrx shape =',block.shape)  # inspect the array shape


# In[15]:


'''Use Matrix Multiplication'''

shape = block.shape[0:2]                  # use the shape to automate noise_mtrx generation
noise_mtrx = np.random.random(shape)      # generate random matrix
block_noise = block[:,:,2] * noise_mtrx   # apply noise to the blue channel
plt.figure(2)
plt.imshow(block_noise)
plt.title('Noisy Matrix',fontsize=14)


# In[16]:


'''Matrix Scaling (matrix product or division by a scalar)'''

mat1 = np.random.random((3,3))
print('mat1     =\n',mat1)
factor = 3.21
result = factor * mat1        # scaling of mat1 element-by-element; product with factor
print('product   =\n',result)


# In[17]:


'''Matrix Scaling of an Image'''

color_channel = np.copy(block[:,:,0])   # copy the red channel
color_channel /= color_channel.max()    # scale to gray, 0-255 values
color_channel *= 255
gray_channel  = color_channel.astype(int) # truncate all float data type to int

plt.figure(3)
plt.imshow(gray_channel,cmap='gray')
plt.title('Matrix Gray Scaling',fontsize=14)


# In[18]:


'''Other Mathematical Operations on a Matrix'''

mtrx = np.copy(block[:,:,0])
plt.figure(4)
plt.imshow(mtrx)
plt.title('Original',fontsize=14)

log_mtrx = np.log(mtrx + .001)  # natural log element-by-element
plt.figure(5)
plt.imshow(log_mtrx)
plt.title('Log Transform',fontsize=14)

exp_mtrx = np.exp(log_mtrx)     # exponential
plt.figure(6)
plt.imshow(exp_mtrx)
plt.title('Exp of Log Transform',fontsize=14)

sin_mtrx = np.sin( mtrx + np.pi/2 )  # sine
plt.figure(7)
plt.imshow(sin_mtrx)
plt.title('Sine Transform',fontsize=14)

mtrx_cubed = mtrx**3            # powers
plt.figure(8)
plt.imshow(mtrx_cubed)
plt.title('Cube Transform',fontsize=14)

mtrx_mean = mtrx.mean()         # arithmetic mean
print('mean(mtrx) =',mtrx_mean)

mtrx_std = mtrx.std()           # standard deviation
print('std(mtrx) =',mtrx_std)


# In[19]:


'''Matrix Transposition'''
'''clockwise rotation followed by horizontal right to left flip'''

mtrx = np.random.random((5,7))
np.set_printoptions(precision=3) # one way to control printing of numpy arrays
print('mtrx =\n',mtrx)

mtrx_T = mtrx.transpose()       # transpose of a mtrx: M[i,j] -> M[j,i]
print('mtrx^T =\n',mtrx_T)


# In[20]:


'''Matrix Transposition'''
'''example of adding a transformed matrix to another transform transposed'''
'''note: to add a matrix to its transpose, a matrix must be square'''

n_rows = block.shape[0]
n_columns = n_rows

mtrx = np.copy(block[:n_rows,:n_columns,0])   # select a square block

sin_mtrx = np.sin( mtrx + np.pi/2 )  # sine
sin_mtrx /= sin_mtrx.max()
plt.figure(9)
plt.imshow(sin_mtrx)
plt.title('Sine Transform',fontsize=14)

mtrx_cubed = mtrx**3                # powers
plt.figure(10)
plt.imshow(mtrx_cubed)
plt.title('Cube Transform',fontsize=14)

plt.figure(11)
plt.imshow( sin_mtrx + mtrx_cubed.transpose() )    # sine + cubed transposed
plt.title('Sine + Cube Transpose Transform',fontsize=14)


# In[ ]: