%Name: Assignment 4 Green Swarm and NE Outlier Eartquakes
%
%This script loads in the Yellowstone Earthquake Swarm of 2009 from a text
%file, plots the data and then saves the file as a jpg titled
%"Yellowstone_EQ".  It then makes figures of the Yellowstone earthquakes
%plotted by lat and long.  There are also options for boxing in an area of earthquakes, naming that area and 
%making a histogram of earthquakes per year in order to find out more
%specific information about the designated earthquakes.
%
%Wendy Bohon, modified from Lecture 4 Computers in Geology class 598
%Last modified 2/18/2009

%Loads in txt file.  Assigns the variables
load('ystone_eqs_sort.txt','-ascii');
year = ystone_eqs_sort(:,1);
x_locations = ystone_eqs_sort(:,4);
y_locations = ystone_eqs_sort(:,3);
depth = ystone_eqs_sort(:,5);
magnitudes = ystone_eqs_sort(:,6);

%This plots the figure, and the clf clears whatever was in the figure
%before.  This also labels the x and y variables and sets the axis as equal
figure(1)
clf
plot(x_locations,y_locations,'k.')
xlabel('longitude')
ylabel('latitude')
axis equal

%this tells the user what to do
fprintf(1, 'click on lower left and upper right of the area of interest\n')

%this tells the user to click on the right and left corners to make a box.
%\n means add a new line.  The argument of 2 tells the ginput "2 clicks"
fprintf(1, 'click on lower left and upper right of the area of interest\n')
[X,Y] = ginput(2);

xmin=min(X)
xmax=max(X)
ymin=min(Y)
ymax=max(Y)

%this draws a box. Realize that it has 5 vertices because you have to close it back around to%this draws a box. Realize that it has 5 vertices because you have to close it back around to the beginning
xs = [xmin xmax xmax xmin xmin];
ys = [ymin ymin ymax ymax ymin];

hold on
plot(xs, ys, 'g-') 

%This is an interactive command that tells you about the events in the box
%tf is for true/false
tf = x_locations <= xmax & x_locations >= xmin & y_locations <= ymax & y_locations >= ymin;
locs = find(tf);
plot(x_locations(locs),y_locations(locs),'ro')

%this puts a label on the box.  Input writes the message for the user.  The
%user can choose the color of the label
your_label = input('What is your label for that selection? \n', 's');

gtext(your_label, 'Color', 'g')

%Makes histograms of event magnitudes.  Labels axes
figure(3)
clf
hist(magnitudes(locs))
xlabel('Magnitude')
ylabel('# events')

%Plots histogram of magnitudes and number of EQ per year(2 windows in one figure)
subplot(2,1,1)
hist(magnitudes(locs))
xlabel('Magnitude')
ylabel('# events')

%This is the second window in the histogram figure 3 plot
subplot(2,1,2)
n = hist(year(locs), years_for_hist);
bar(years_for_hist, n)
xlabel('Event year')
ylabel('# events')

%This script load in the Yellowstone Earthquake Swarm of 2009 from a text
%file, plots the data and then saves the file as a jpg titled
%"Yellowstone_EQ".  This script is the same as above, but modified to
%select a different group of earthquakes.

%Loads in txt file
load('ystone_eqs_sort.txt','-ascii');
year = ystone_eqs_sort(:,1);
x_locations = ystone_eqs_sort(:,4);
y_locations = ystone_eqs_sort(:,3);
depth = ystone_eqs_sort(:,5);
magnitudes = ystone_eqs_sort(:,6);

%This plots the figure, and the clf clears whatever was in the figure
%before
figure(4)
clf
plot(x_locations,y_locations,'k.')
xlabel('longitude')
ylabel('latitude')
axis equal

fprintf(1, 'click on lower left and upper right of the area of interest\n')

%this tells the user to click on the right and left corners to make a box.
%\n means add a new line
fprintf(1, 'click on lower left and upper right of the area of interest\n')
[X,Y] = ginput(2);

xmin=min(X)
xmax=max(X)
ymin=min(Y)
ymax=max(Y)

%this draws a box. Realize that it has 5 vertices because you have to close it back around to%this draws a box. Realize that it has 5 vertices because you have to close it back around to the beginning
xs = [xmin xmax xmax xmin xmin];
ys = [ymin ymin ymax ymax ymin];

hold on
plot(xs, ys, 'm-') 

%This is an interactive command that tells you about the events in the box
tf = x_locations <= xmax & x_locations >= xmin & y_locations <= ymax & y_locations >= ymin;
locs = find(tf);
plot(x_locations(locs),y_locations(locs),'mo')

%this puts a label on the box.  Input writes the message for the user
your_label = input('What is your label for that selection? \n', 's');

gtext(your_label, 'Color', 'm')

%Makes histograms of event magnitudes
figure(5)
clf
hist(magnitudes(locs))
xlabel('Magnitude')
ylabel('# events')

%Plots histogram of magnitudes and number of EQ per year(2 windows in one figure)
subplot(2,1,1)
hist(magnitudes(locs))
xlabel('Magnitude')
ylabel('# events')

subplot(2,1,2)
n = hist(year(locs), years_for_hist);
bar(years_for_hist, n)
xlabel('Event year')
ylabel('# events')