clear all
close all
figure();
hold on;
format longEng



lat_min = 61; %center is (65,-19)
lat_max = 69;
lon_min = -34;
lon_max = -4;
axis([lon_min lon_max lat_min lat_max]);



ocean = [lat_min, lon_min; lat_min, lon_max; lat_max, lon_max; lat_max, lon_min; ];
color1 = [0 0.3 1];
h = patch(ocean(:,2),ocean(:,1),color1);



faroeislands = load('faroeislands.txt');
color2 = [0 0 1];
h = patch(faroeislands(:,2),faroeislands(:,1),color2);

greenland = load('greenland.txt');
color3 = [0 0 1];
h = patch(greenland(:,2),greenland(:,1),color3);

iceland = load('iceland.txt');



v1 = [66.3000, -17.1000; 
      65.8944, -16.8406;
      65.7171, -16.7544;
      65.0111, -16.7485;
      64.6530, -16.6470;
      64.6410, -17.5280; 
      64.4167, -17.3222;
      64.18054, -18.24386;
      64.0707, -18.2378;
      63.9664, -18.6082;
      63.9923, -19.6658;
      64.0686, -20.6431;
      64.0769, -21.1997;
      63.8919, -22.0665;
      63.8900, -22.5200;
      63.6700, -23.3300;]; %coordinates of volcanos in mid-atlantic ridge

  
  
WesternIceland = load('leftsideiceland.txt');
EasternIceland = load('rightsideiceland.txt');



v2 = [64.1265, -21.8174; 
      65.6885, -18.1262;
      63.9998, -22.5583;
      64.3219, -22.0703;
      66.0611, -23.1889;
      65.2669, -14.3948;
      63.9335, -21.0014;];  %coordinate of most populate cities in Iceland per region
 
 p2 = size(v2);
 j2 = p2(1);
 
  txt ={'\leftarrow Reykjavík', 
        '\leftarrow Akureyri',
        '\leftarrow Keflavík',
        '\leftarrow Akranes',
        '\leftarrow Ísafjörður',
        '\leftarrow Egilsstaðir',
        '\leftarrow Selfoss'};
  
    
    
%rate of 18.2 mm/year and a direction of 105
v = 0.018;
vyleft = v*sind(165);
vxleft = v*cosd(165);
vyright = v*sind(-15);
vxright = v*cosd(-15);

nyleft = size(WesternIceland,1);
nxleft = size(WesternIceland,2);
nyright = size(EasternIceland,1);
nxright = size(EasternIceland,2);

R = 6371000;
latt = 2*pi*R/360;
long = pi*R/180*cosd(65);

WesternIcelandNew = zeros(nyleft,nxleft);
EasternIcelandNew = zeros(nyright,nxright);
v2new = zeros(7,2);

t = 0; %Change years from 0 to 25000000 in increments of 1000000

for i = 1:nyleft
aleft = [WesternIceland(i,1),WesternIceland(i,2)]+[(vyleft*t)/latt,(vxleft*t)/long];
WesternIcelandNew(i,:) = WesternIcelandNew(i,:) + aleft;
end
color4 = [0.4 0.2 0.2];
h = patch(WesternIcelandNew(:,2),WesternIcelandNew(:,1),color4);

for i = 1:nyright
aright = [EasternIceland(i,1),EasternIceland(i,2)]+[(vyright*t)/latt,(vxright*t)/long];
EasternIcelandNew(i,:) = EasternIcelandNew(i,:) + aright;
end
color5 = [0.8 0.4 0.4];
h = patch(EasternIcelandNew(:,2),EasternIcelandNew(:,1),color5);
 
 for i = 1:5
 b1 = [v2(i,1), v2(i,2)] + [(vyleft*t)/latt,(vxleft*t)/long];
 v2new(i,:) = v2new(i,:) + b1;
 end
 
 for i = 6:7
  b2 = [v2(i,1), v2(i,2)] + [(vyright*t)/latt,(vxright*t)/long];
  v2new(i,:) = v2new(i,:) + b2;
 end

for count = 1:5
    plot(v2new(count,2),v2new(count,1),'og');
    text(v2new(count,2),v2new(count,1),txt(count));
end

for count = 6:7
    plot(v2new(count,2),v2new(count,1),'og');
    text(v2new(count,2),v2new(count,1),txt(count));
end

num = num2str(t);
title(['t = ', num, ' years']) ;