'''
vander_field field ploter for vanderpol driven oscillator
writen by Timothy Jones, Drexel University, 2010

Copyright (C) 2010 Timothy Jones

    This program is free software: you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    This program is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with this program.  If not, see <http://www.gnu.org/licenses/>.


'''
from __future__ import division
from pylab import *
from numpy import ma
import numpy

pi =  3.1415926535
a  =  0.25
b  =  0.7
c  =  1.0
d  = 10.0
omega = pi/2

def func1(x,y):
    return  b*y+(c-d*y*y)*x
def func2(x,y,t):
    return  -x + a*sin(omega*t)

ii=0.0
counter=0
while(ii<=4.0):

    
    num1=ii*1.0
    num2=100000+counter

    file1="./slices/%.6f" % num1
    figname="./images/%.6f.png" % num1

    file2='./slices2/%i.dat' % num2

    imported_array = numpy.loadtxt(file1,delimiter=' ')
    #print imported_array.shape
    xx1 = imported_array[:,0]
    yy1 = imported_array[:,1]
    imported_array2 = numpy.loadtxt(file2,delimiter=' ')
    #print imported_array2.shape
    xx2 = imported_array2[:,0]
    yy2 = imported_array2[:,1]

   

  

    X,Y = meshgrid( arange(-1,1,.08),arange(-1,1,.08) )
    XX,YY = meshgrid( arange(-1,1,.02),arange(-1,1,.02) )

    A = func1(XX,YY);
    B = func2(XX,YY,ii);
    C = func1(X,Y);
    D = func2(X,Y,ii);


    figure(1)
  
    plot(xx1,yy1,'o',markerfacecolor='white', markeredgecolor='red')
    plot(xx2,yy2,'.',markersize=2)
    M = log(0.2+sqrt(pow(A, 2) + pow(B, 2)))


    '''
    pcolor(XX,YY,M)
    colorbar([-6,2])
    '''

    #ax = subplot(111)
    im = imshow(M, cmap=cm.jet,extent=(-1,1, 1,-1))
    #im.set_interpolation('nearest')
    #im.set_interpolation('bicubic')
    im.set_interpolation('bilinear')
    #colorbar()

    Q = quiver( X, Y, C, D, units='x') 
    

    axis([-1, 1, -1, 1])
    savefig(figname)
    close(figure(1))
    ii+=0.02
    counter+=1
    #show()