program analysis;
{
This program calculates all available station derivatives each
year and their average to estimate the global derivative. First
the program reads in the lines of GYA.txt, the yearly averages for
all station. It sorts these in order of station and year. It turns
out that our data file was pre-sorted.
}
const
max_num_measurements=1000000;
type
measurement_type=record
id:longint;{station identifier}
year:integer;{year of the measurement}
temperature:real;{year-long average}
derivative_valid:boolean;{valid if last year's measurement exists}
derivative:real;{change from last year}
end;
var
data:text;{the data file}
measurements:array [1..max_num_measurements] of measurement_type;
measurement_num,num_measurements,y,counter:integer;
i:longint;
t,sum:real;
sorted:boolean;
{
We use this swap procedure to sort our list of stations. We
use a "bubble sort" because it keeps the code simple.
}
procedure swap(p,q:integer);
var
m:measurement_type;
begin
m:=measurements[p];
measurements[p]:=measurements[q];
measurements[q]:=m;
end;
{
This is the main program.
}
begin
{
Start by reading in the Global Yearly Averages from our GYA.txt file,
which we have formatted for our convenience, but which contains raw
data from NCDC. Each line in the file provides us with a station ID,
a year, and a temperature.
}
writeln('Read GYA.txt...');
reset(data,'GYA.txt');
num_measurements:=0;
while not eof(data) do begin
inc(num_measurements);
with measurements[num_measurements] do
readln(data,id,year,temperature);
end;
close(data);
{
We sort the measurements by increasing station ID and increasing
year.
}
writeln('Sort measurements by station and year...');
sorted:=false;
while not sorted do begin
sorted:=true;
for measurement_num:=1 to num_measurements-1 do begin
if (measurements[measurement_num].id > measurements[measurement_num+1].id)
or
((measurements[measurement_num].id = measurements[measurement_num+1].id)
and
(measurements[measurement_num].year > measurements[measurement_num+1].year) )
then begin
swap(measurement_num,measurement_num+1);
sorted:=false;
end;
end;
end;
{
Each measurement now gets a derivative, provided that the previous measurement
is from the same station in the previous year.
}
writeln('Calculate and mark all valid derivatives...');
measurements[1].derivative_valid:=false;
for measurement_num:=2 to num_measurements do begin
if ((measurements[measurement_num].id = measurements[measurement_num-1].id)
and
(measurements[measurement_num].year - 1 = measurements[measurement_num-1].year) )
then begin
measurements[measurement_num].derivative_valid:=true;
measurements[measurement_num].derivative:=
measurements[measurement_num].temperature
- measurements[measurement_num-1].temperature;
end else begin
measurements[measurement_num].derivative_valid:=false;
measurements[measurement_num].derivative:=-999;
end;
{ with measurements[measurement_num] do writeln(id:1,' ',year:1,' ',derivative:10:6);}
end;
{
For each year in our time range, calculate the average derivative of temperature
and print it to the screen.
}
writeln('Calculate average derivative for each year...');
for y:=1800 to 2010 do begin
sum:=0;
counter:=0;
for measurement_num:=1 to num_measurements do
with measurements[measurement_num] do
if (year=y) and derivative_valid then begin
sum:=sum+derivative;
inc(counter);
end;
writeln(y:1,' ',sum/counter:10:6,' ',counter:1);
end;
{
What you see now in the console is the derivatives with year. We pasted these into
Excel and integrated them. You will find our Climate.xls spreadsheet that performs the
integration in our Climate.zip archive.
http://www.hashemifamily.com/Kevan/Climate/Climate.zip
As you can see, the integration proceeds by adding the previous year's integral to
the current year's average derivative. We perform no smoothing.
}
end.