Saturday, July 31, 2010

The story of Nachiketa

This has been picked up from Dr K S Narayanachar's book..And i must say this story has really started giving me a sense of satisfaction..i will be sharing some of it here..
This story appears primarily in Kathopanishad..It is also present in RigVeda 10th mandala..Also in
Taittiriya Brahmana, in abridged form..
Nachiketa was the son of a rishi named Vaajashrava..Vaja means food..Shrava means fame..Vaajashrava
had earned a lot of fame through "Annadaan"(donating food or feeding people)..The beauty of our hindu scriptures is that names of characters themselves tell a lot about their culture..Now after some time this fame probably started having its side effects..vaajashrava decided to perform a great
yagna called "vishwajith"..A performing this should ideally give away everything he possesses(including himself)..But vaajashrava was probably after "fame" than god..In fact the main reason for performing yagnas is to develop the habit of donating clothes,food etc and to make one understand that "nothing in this world is mine"..Vaajashrava was giving away old, useless cows as donation to people..This child nachiketa was watching his father..He realized that the yagna was becoming a failure..He went to his father and said "To whom will u donate me??"..The sage simply ignored this child, but nachiketa continued asking again and again. Vaajashrava lost his patience and said "I ve given u as daan to yama,the god of death".

Here the rig veda and also kathopanishad describe in a beautiful way how nachiketa felt and what did he think..Unlike us, nachiketa did not get angry at this..He thought "My father is a great sage..whatever he says is for my good!! But of what use am i to yama??..even plants eat,grow,reproduce and perish..I think i am not the worst among his disciples..So now that i ve been given as offering to yama, of what use i can be to him??"..of course i cant fully describe the opinion of vedas as in the book..i dont think english is suitable for this..
But the point here is that the child was able to think positively..Our indian culture teaches us to respect elders..even if they go wrong, youngsters should not speak in a manner which hurts their sentiments..nachiketa was a firm believer in this..This belief took him to the abode of god..

Thursday, July 29, 2010

video shot detection and summary

This is my first python project..not upto coding standards ofcourse..
#######################################################################################

import cv
import numpy
import time
from threading import Thread
from Queue import Queue
#now that we ve failed with hsv let us try rgb

hist_size = 256 #the number of bins
range_0 = [0, 256] #ranges of values to be stored in bins(here 8 bits)
ranges = [ range_0 ]

class Worker(Thread):

def __init__(self):
Thread.__init__(self)

self.queue = Queue()
self.writer = None
# Daemon threads won't prevent process from exiting
self.setDaemon(True)

# Start ourselves automatically
self.start()

def run(self):
writer = None
while 1:
frame = self.queue.get()

if not frame:
# Asked to exit
return

# Create a writer if needed - implement as necessary.
# (add init parameters, methods to change/reset, etc...)
# This assumes the frame size won't change
if not writer:
writer = cv.CreateVideoWriter("output.mp4", cv.CV_FOURCC('D','I','V','X'), 25, (640,480),1)

cv.WriteFrame(writer, frame)

# Requests from main thread
def write(self, frame):
# Duplicate to ensure we hang onto an independent image
self.queue.put(cv.CloneImage(frame))

def stop(self):
self.queue.put(None)

def rgb_hist(src_img):
#one histogram for each channel
hist_red = cv.CreateHist([hist_size], cv.CV_HIST_ARRAY, ranges,1)
hist_blue = cv.CreateHist([hist_size], cv.CV_HIST_ARRAY, ranges,1)
hist_green= cv.CreateHist([hist_size], cv.CV_HIST_ARRAY, ranges,1)
blue_img = cv.CreateImage( cv.GetSize(src_img), cv.IPL_DEPTH_8U,1)
green_img = cv.CreateImage( cv.GetSize(src_img),cv.IPL_DEPTH_8U,1)
red_img = cv.CreateImage(cv.GetSize(src_img),cv.IPL_DEPTH_8U,1)
cv.Split( src_img, blue_img, green_img, red_img,None)
#now calculate histogram for each channel
cv.CalcHist([cv.GetImage(red_img)], hist_red)
cv.CalcHist([cv.GetImage(blue_img)], hist_blue)
cv.CalcHist([cv.GetImage(green_img)], hist_green)
##allright!! finally i am able to read the rgb stuff
##find_max = cv.QueryHistValue_1D(hist_red,200)
##print find_max
return hist_red,hist_blue,hist_green


def calcHistDiff(hist1,hist2):
sum=0
NZ1=0
NZ2=0
red_hist1=hist1[0]
blue_hist1=hist1[1]
green_hist1=hist1[2]

red_hist2=hist2[0]
blue_hist2=hist2[1]
green_hist2=hist2[2]
for i in range(256):
red_binvalue1=cv.QueryHistValue_1D(red_hist1,i)#get bin values
red_binvalue2=cv.QueryHistValue_1D(red_hist2,i)#get bin values
if red_binvalue1:
NZ1 = NZ1 + 1

if red_binvalue2:
NZ2 = NZ2 + 1


blue_binvalue1=cv.QueryHistValue_1D(blue_hist1,i)#get bin values
blue_binvalue2=cv.QueryHistValue_1D(blue_hist2,i)#get bin values
if blue_binvalue1:
NZ1 = NZ1 + 1

if blue_binvalue2:
NZ2 = NZ2 + 1


green_binvalue1=cv.QueryHistValue_1D(green_hist1,i)#get bin values
green_binvalue2=cv.QueryHistValue_1D(green_hist2,i)#get bin values
if green_binvalue1:
NZ1 = NZ1 + 1

if green_binvalue2:
NZ2 = NZ2 + 1

#now increment respective non-zero bin counts
#if red_binvalue1 or blue_binvalue1 or green_binvalue1:
# NZ1 = NZ1 + 1

#if red_binvalue2 or blue_binvalue2 or green_binvalue2:
#NZ2 = NZ2 + 1

sum=sum+abs((red_binvalue1 + blue_binvalue1 + green_binvalue1) - (red_binvalue2 + blue_binvalue2 + green_binvalue2))

ADH=sum/256.0
return ADH,NZ1,NZ2

def calcWeight(NZ,window):
NZM=0.0
weight=0.0
nz_sum=numpy.sum(NZ[:-1])
#assert nz_sum > 0,'this cant be zero'
#if nz_sum==0:
# nz_sum=1
NZM=(nz_sum / window)
weight=(512.0)/NZM

#cv.WaitKey(0)
return weight*3

def calcThreshold(ADH,curr_frame,NZ):
threshold=0.0
weight=0.0
adh_sum=numpy.sum(ADH[:-1])#since first value is ADH of previous shot frame
#window=curr_frame - prev_shot_frame
window=len(ADH)-1

weight=calcWeight(NZ,window)
#weight=2.7
threshold=(adh_sum/(window + 1))*weight
#print adh_sum
#print threshold
#cv.WaitKey(0)
return threshold

def calcOptiFlow(img1,img2):
#first convert to grayscale
gray_img1=cv.CreateImage((img1.width,img1.height),cv.IPL_DEPTH_8U,1)
gray_img2=cv.CreateImage((img1.width,img1.height),cv.IPL_DEPTH_8U,1)
cv.CvtColor(img1,gray_img1,cv.CV_RGB2GRAY)
cv.CvtColor(img2,gray_img2,cv.CV_RGB2GRAY)

#create velocity matrices for both x and y direction
dx = cv.CreateMat(img1.height,img1.width,cv.CV_32FC1)
dy = cv.CreateMat(img1.height,img1.width,cv.CV_32FC1)

#now calculate opticalflow
cv.CalcOpticalFlowLK(gray_img1,gray_img2,(5,5),dx,dy)

#now calculate bulk velocity in both directions
bulk_x = cv.Avg(dx)
bulk_y = cv.Avg(dy)
return bulk_x[0],bulk_y[0]


def calcOpticFlow(img1,img2):
max_count=500
min_distance=0.01
quality=0.01
win_size=3
dx=0
dy=0
flags=0
features_img1=[]
features_img2=[]
status=[]
track_error=[]
MAX_COUNT=200
motion_intensity=0
orientation=0
#first convert to grayscale
gray_img1=cv.CreateImage((img1.width,img1.height),cv.IPL_DEPTH_8U,1)
gray_img2=cv.CreateImage((img2.width,img2.height),cv.IPL_DEPTH_8U,1)
cv.CvtColor(img1,gray_img1,cv.CV_RGB2GRAY)
cv.CvtColor(img2,gray_img2,cv.CV_RGB2GRAY)
#allocate other stuff
eig = cv.CreateImage (cv.GetSize (gray_img1), 32, 1)
temp = cv.CreateImage (cv.GetSize (gray_img2), 32, 1)

#LK carves images imto pyramids
pyramid1 = cv.CreateImage (cv.GetSize (img1), 8, 1)
pyramid2 = cv.CreateImage (cv.GetSize (img2), 8, 1)


#get good features of both images
features_img1 = cv.GoodFeaturesToTrack (gray_img1, eig, temp, MAX_COUNT, quality, min_distance, None, 3, 0, 0.04)
#features_img2 = cv.GoodFeaturesToTrack (gray_img2, eig, temp, MAX_COUNT, quality, min_distance, None, 3, 0, 0.04)
features_img2,status,track_error=cv.CalcOpticalFlowPyrLK(gray_img1,gray_img2,pyramid1,pyramid2,features_img1,(win_size, win_size), 3,(cv.CV_TERMCRIT_ITER|cv.CV_TERMCRIT_EPS, 20, 0.03), flags)
features_img2 = [ p for (st,p) in zip(status, features_img2) if st]
#print features_img2
#cv.WaitKey(0)
for i in range(len(features_img2)):
x1=features_img1[i][0]
y1=features_img1[i][1]
x2=features_img2[i][0]

y2=features_img2[i][1]

dx=dx+abs(x1-x2)
dy=dy+abs(y1-y2)

dx=dx/4
dy=dy/4
motion_intensity=motion_intensity+numpy.sqrt(numpy.square(dx) + numpy.square(dy))
if dx > 0.0:
orientation=orientation+numpy.arctan(dy/dx)
return (motion_intensity,orientation)

#this is the final and most important step ie creating video file
def writeVideo(keyframe_list):
i=0
input_video = "test_videos/manipal.MP4"
len_arr=len(keyframe_list)
worker=Worker()
cap=cv.CreateFileCapture(input_video)
img1=cv.QueryFrame(cap)
while 1:
i=i+1
curr_frame=i+1
img2=cv.QueryFrame(cap)
if not img2:
break
if not len_arr:
print "dude no keyframes"
break
if curr_frame > keyframe_list[0]:
#we are here means that we have already written the keyframe
keyframe_list=keyframe_list[1:]
len_arr=len_arr - 1
continue

if curr_frame == keyframe_list[0]:
worker.write(img2)
#also write 3 extra frames
img2=cv.QueryFrame(cap)
worker.write(img2)

img2=cv.QueryFrame(cap)
worker.write(img2)

img2=cv.QueryFrame(cap)
worker.write(img2)
keyframe_list=keyframe_list[1:]
len_arr = len_arr -1
i=i+2

k = cv.WaitKey(5)
if k % 0x100 == 27:
# ESC key to exit
break
# Shut down worker thread and wait for it to exit in case it's
# in the middle of doing some I/O
worker.stop()
worker.join(5)
if __name__=='__main__':
#img=cv.LoadImage("/home/guruvyasa/openCV_samples/lena.jpg")
#rgb_hist(img)
input_video = "test_videos/manipal.MP4"
prev_shot_frame=0#assuming first frame as shot change frame
curr_frame=0
shot_frame_list=[]#list holds framenos of shot frames
val=0
ADH=[0]
temp=0
#x=1#keeps track of index in NZ list
NZ=[0]
NZ1=0
NZ2=0
motion_intensity=[]
keyframe_list=[]
#cv.NamedWindow("first",1)
cap=cv.CreateFileCapture(input_video)
img1=cv.QueryFrame(cap)
#image1=cv.LoadImage("/home/guruvyasa/openCV_samples/lena.jpg")
#image2=cv.LoadImage("/home/guruvyasa/openCV_samples/baboon.jpg")
#hist1=rgb_hist(image1)
#hist2=rgb_hist(image2)
#ADH=calcHistDiff(hist1,hist2)
hist1=rgb_hist(img1)
i=0

for i in range(700):#the frame no is i+1
#while 1:
i=i+1

curr_frame=i+1
#if (curr_frame == 60) or (curr_frame == 165) or (curr_frame ==220):
# print "********i am here*******"
# cv.WaitKey(0)
img2=cv.QueryFrame(cap)
if not img2:
break
#hist1=rgb_hist(img1)
hist2=rgb_hist(img2)
if ADH > 0:#this is always true..i am lazy not to delete this one
temp,NZ1,NZ2=calcHistDiff(hist1,hist2)
ADH.append(temp)
#testing
#print x
NZ[-1]=NZ1
NZ.append(NZ2)
#now calcuate the adaptive threshold
#print ADH
#print NZ
#if curr_frame > 270 and curr_frame < 293:
# print motion_intensity
# cv.WaitKey(0)
threshold=calcThreshold(ADH,curr_frame,NZ)
#print threshold
cv.SaveImage("temp1.jpg",img1)
cv.SaveImage("temp2.jpg",img2)
img1=cv.LoadImage("temp1.jpg")
img2=cv.LoadImage("temp2.jpg")
block=((img1.width)/2,(img1.height)/2,(img1.width)/2,(img1.height)/2)
blk_img1=cv.GetSubRect(img1,block)
blk_img2=cv.GetSubRect(img2,block)

val1,val2=calcOpticFlow(blk_img1,blk_img2)
motion_intensity.append(val1)
#if it has been quite some time since a shot(hardcoded)
if len(ADH) == 60:
#now go on to calculate maximum motion intensity in the shot
m_max=numpy.max(motion_intensity)
#print m_max
#now we take only those frames with high motion intensity
temp_x=curr_frame-2
len_adh=len(ADH)-5 #dont consider current frame
for z in range(2,len_adh):
#temp_x=temp_x-1
temp_val1=motion_intensity[-z] #traverse the list in reverse
temp_val2=motion_intensity[-(z+1)]
if abs(temp_val2-temp_val1) > 1500:
#print temp_x
#cv.WaitKey(0)
keyframe_list.append(temp_x)
temp_x=temp_x-1
#now since we have processed ADH for keyframes we resize it
ADH=ADH[-1:]
NZ=NZ[-1:]

if (ADH[-1] > threshold) and ADH[-1] > 10.0:
shot_frame_list.append(curr_frame)
keyframe_list.append(curr_frame)
#key_frames.append(curr_frame)
cv.ShowImage("prev",img1)
cv.ShowImage("current",img2)
cv.WaitKey(0)
ADH=ADH[-1:]#here we clear the list till current frame
#so now the first value in ADH and NZ is value of previus shot frame
NZ=NZ[-1:]
#x=1

img1=img2
hist1=hist2

#print len(keyframe_list)
#cv.WaitKey(0)
keyframe_list=numpy.sort(keyframe_list)
#print keyframe_list
#now its time for showdown, write the frames for summary
#writeVideo(keyframe_list)

Eigen values and eigen vectors

i always had this doubt as to what these eigen stuffs are..till i went through this http://ceee.rice.edu/Books/LA/eigen/

Calculus continued

now, in differentiation we need to understand what dx,dy etc mean actually..they are called differentials..ie infinitesimally small increments in x,y etc..while programming we can make this as a variable say step_size..
this simple knowledge leads us into integration..
Integration is simply inverse of differentiation..when we integrate we are calculating areas under the curve..
Why do we need to do that?? well, suppose we know that dv(infinitesimally small change in velocity) is related to dt(time change) as dv=f(t)dt ie there is a function which tells us how a small change in time affects velocity we can calculate the relation between velocity and time over an interval as a whole using integration..
now, what does f(t)dt mean?? it s nothing but the area of a rectangle with infinitesimally small width..the smaller the width, better is the approximation for the area

Calculus stuff

I ve started studying calculus again..The last time i did this, there were many questions and i backed out..Now i am able to understand many things and hence i ll share some things here..
1) Differetiation is all about slopes..Slopes of tangents to curves..Now i had always wondered as to why we need to worry about curves and their slopes..After studying many topics with these, now i realize that we can model almost anything using curves(equations)..ie mathematics is like a programming language. We need to represent real-world problems using some notations..eg..suppose we know(from experiments and observations) that pixel intensity(in images) depends on two factors(say x and y)..we can assume a linear dependence and say pixel_intensity=m*x + n*y..thus as u can see we have an equation of a line now..
now, slope=(change in y)/(change in x)..ie if we get the slope we can get a relation, of how 'y' changes w.r.t 'x'
Entire differentiation circles around this..it s really amazing..