Dear All,
I have tried to use the code kindly provided by Robert Picard in this post (#10) for precipitation data, to obtain information about air temperature from 1990 to 2014. I used the data available here taken from Willmott, Matsuura and Collaborators' Global CLimate Resource Page and the shapeffile available here. After the construction I kept only year 1900.
However, the final data look unrealistic. The summary statistics are reported below:
Initially, I though that the problem could be in the command highlighted in red above (so instead of summing up through rows, I took the mean - egen p = mean(p1-p12). But also this attempt was not successful (in this case I obtain negative temperature for all countries). Could you please get me suggestions about how to fix the code above?
Thanks for any help you may provide.
Dario
I have tried to use the code kindly provided by Robert Picard in this post (#10) for precipitation data, to obtain information about air temperature from 1990 to 2014. I used the data available here taken from Willmott, Matsuura and Collaborators' Global CLimate Resource Page and the shapeffile available here. After the construction I kept only year 1900.
Code:
shp2dta using "Original data\TM_WORLD_BORDERS-0.3.shp", ///
data("Original data/world_data.dta") coor("Original data/world_coor.dta") replace
use "Original data\world_data.dta", clear
list
save "Original data/allworld_data.dta", replace
use "Original data/world_coor.dta", clear
gen long obs = _n
merge m:1 _ID using "Original data/allworld_data.dta", keep(match) keepusing(_ID) nogen
sort obs
save "Original data/allworld_coor.dta", replace
* Importing weather data one time period
infile lon lat p1 p2 p3 p4 p5 p6 p7 p8 p9 p10 p11 p12 ///
using "Original data\air_temp_2014/air_temp.1900", clear
* Reducing to grid coordinates
keep lat lon
order lat lon
* Matching grid coordinates to countries
geoinpoly lat lon using "Original data/allworld_coor.dta"
* Dropping grid points that did not match the selected countries
drop if mi(_ID)
* Showing the matching grid points
scatter lat lon if _ID == 65 , msize(tiny) mcolor(blue) ///
|| ///
scatter lat lon if _ID != 65 , msize(tiny) mcolor(blue)
graph export "Original data/word2.png", width(800) replace
save "Original data/allworld_grid.dta", replace
clear all
* Obtaining a list of the files in the "air_temp_2014" subdirectory
filelist , dir("Original data\air_temp_2014")
* Extracting the year from the file name
gen year = real(subinstr(filename,"air_temp.","",1))
assert !mi(year)
sum year
* Program to create a file with all the data for each year
program input_each_year
// copy file path and year to local macros
local filename = dirname + "/" + filename
local thisyear = year
dis "---------------------------- processing `filename'"
infile lon lat p1 p2 p3 p4 p5 p6 p7 p8 p9 p10 p11 p12 ///
using "`filename'", clear
gen int year = `thisyear'
egen p = rowtotal(p1-p12)
drop p1-p12
// Merge each grid point to get the _ID of the matched country
merge 1:1 lat lon using "Original data/allworld_grid.dta", ///
assert(master match) keep(match) nogen
// Convert to a long layout
order lat lon year p
isid lat lon year, sort
// What's left at this point is considered results and accumulates
end
runby input_each_year, by(year) status
isid lat lon year, sort
save "Original data/air_temp_annual_grid_1900_to_2014.dta", replace
clear all
infile lon lat p1 p2 p3 p4 p5 p6 p7 p8 p9 p10 p11 p12 ///
using "Original data\air_temp_2014/air_temp.1900", clear
* Reducing to grid coordinates
keep lat lon
order lat lon
isid lat lon, sort
* Decimal degree to radians constant
scalar d2r = c(pi) / 180
* Distance between .5 degrees of latitude, constant longitude
geodist 0 0 .5 0, sphere
scalar ydist = r(distance)
* Distance between .5 degrees of longitude, constant latitude
gen xdist = 0.5 * d2r * cos(lat * d2r) * 6371
* Double-checking using shortest path distance
gen lon1 = lon - .25
gen lon2 = lon + .25
geodist lat lon1 lat lon2, sphere gen(dcheck)
gen diff = abs(xdist-dcheck)
sum diff
* The area of the grid point's Voronoi diagram
gen grid_area = xdist * ydist
sum grid_area
keep lat lon grid_area
save "Original data/grid_area.dta", replace
clear all
use "Original data/allworld_grid.dta", clear
merge 1:1 lat lon using "Original data/grid_area.dta", keep(match) nogen
merge 1:m lat lon using "Original data/air_temp_annual_grid_1900_to_2014.dta", ///
keep(match) nogen
* Using the grid area to weight the air temperature data
gen p_area = p * grid_area
collapse (count) pN=p (sum) p_area grid_area, by(year _ID)
* Mean of air temperature in C degrees
gen p_mean = p_area / grid_area
merge m:1 _ID using "Original data/allworld_data.dta", ///
keepusing(NAME ISO3) assert(match using) keep(match) nogen
rename (NAME ISO3) (country country3)
isid country year, sort
rename p_mean air_temp
label variable air_temp "mean of air temperature in C degrees"
drop pN p_area grid_area
keep if year==1900
save "Original data/air_temp.dta", replace
Code:
Variable | Obs Mean Std. dev. Min Max -------------+--------------------------------------------------------- air_temp | 189 212.0034 117.7725 -422.5823 343.0486
Thanks for any help you may provide.
Dario
Comment