* ifort -132 -assume byterecl get_area_size_grid-1x1.f

	program gridareax

	real theta(181), phi(360), gridarea(180), AREA(180)
	integer i, j, lon, lat
	
	parameter (a = 6.371e+6, a2= a*a, pi0 = 3.14159265/180.) !earth radius = 6.37 x 10^6 m
	
	j = 0
	i = 0

	do lat = -90, 90, 1
	j = j + 1
	theta(j) = pi0*lat
	write(6,*)'lat j = ', lat, j, theta(j)
	enddo

	do lon =   0, 359, 1
	i = i + 1
	phi(i)  = pi0*lon
	write(6,*)'lon i = ', lon, i
	enddo

	sum = 0	
	delta_phi = phi(2) - phi(1) 
	do j = 1, 180
c	gridarea(j) = a2*(sin (theta(j+1) )- sin( theta(j)))*(phi(i+1) - phi(i))
	gridarea(j) = a2*(sin (theta(j+1) )- sin( theta(j)))*(delta_phi)
	sum = sum + gridarea(j)
	enddo
c	do j = 11, 170	!1 TO 80
c	area(j-10) = gridarea(j)	
c	write(6,*)'j delta_phi theta(j) gridarea(j) = ',j-10, delta_phi, theta(j), gridarea(j)
c	enddo

	do j = 1, 180	!1 TO 80
	area(j) = gridarea(j)	
	write(7,'(a,i3,1x,i3,3e15.4)')'j delta_phi theta(j) gridarea(j) = ',j-91,j-90, delta_phi, theta(j), gridarea(j)
	enddo

	open(unit=10, file='grid_box_area_1x90_1x1.datt', access = 'direct',
     @	form='unformatted', status='unknown', recl = 4*1*180)
	open(unit=50, file='grid_box_area_1x90_1x1.asc',
     @	form='formatted', status='unknown')

	write(10,rec=1)area
	CLOSE (10)
	
	do j = 180, 1, -1 
	jj = j-90
	if (j .LE. 90) jj = -(90 - j + 1)
	
	write(50, '(2i5, 1pe14.6)') j, jj, area (j)
	enddo
	
	sum = sum * 360
	write(6, *)'sum true = ', sum, 4.*3.14159265*a2

	end
	
