pkg/ecco/cost_geoid.F

C $Header: /u/gcmpack/MITgcm/pkg/cost/Attic/cost_geoid.F,v 1.1.2.2 2002/04/04 10:58:58 heimbach Exp $

#include "COST_CPPOPTIONS.h"

      subroutine cost_geoid(
     O                     sphcost,
     I                     shc,
     I                     mythid
     &                   )

c     ==================================================================
c     SUBROUTINE cost_geoid
c     ==================================================================
c
c     o Evaluate the cost function of the geoid contribution.
c
c     started: Christian Eckert eckert@mit.edu 30-Jun-1999
c
c     changed: Christian Eckert eckert@mit.edu 25-Feb-2000
c
c              - Restructured the code in order to create a package
c                for the MITgcmUV.
c
c     changed: Ralf Giering Ralf.Giering@FastOpt.de 12-Jun-2001
c
c              - totally rewrite for parallel processing
c
c     ==================================================================
c     SUBROUTINE cost_geoid
c     ==================================================================

      implicit none

c     == global variables ==

#include "EEPARAMS.h"
#include "SIZE.h"
#include "PARAMS.h"

#include "ecco_cost.h"
#include "cost_sph.h"

c     == routine arguments ==

      _RL     sphcost
      _RL     shc( ncShc )
      integer mythid

#ifdef ALLOW_EGM96_ERROR_COV

c     == local variables ==

      Real*8  pinv( ncshc )

      integer i,j
      integer ilo,ihi
      integer ilnblnk
      integer ifnblnk
      integer ireads
      integer egmunit
      integer egmsize
      parameter( egmsize = 8 )

      integer jsize, joff, jbeg, jend, j

      _RL factor
      _RL recip_rr

c     == external functions ==
      integer  ilnblnk
      external ilnblnk
      integer  ifnblnk
      external ifnblnk

c     == end of interface ==

c--   Only the master thread is doing I/O
      _BEGIN_MASTER( mythid )

c--   initialise cost variable for all threads
c--   to allow usage of global_sum_r8
      sphcost = 0. _d 0
      do i=1,ncShc
         pinv(i) = 0. _d 0
      enddo

c--   Initialise variables.
      recip_rr = recip_rsphere*recip_rsphere

c--   Open the geoid error covariance data file.
      call mdsfindunit( egmunit, mythid )
      ilo = ifnblnk(geoid_covariancefile)
      ihi = ilnblnk(geoid_covariancefile)

      print*
      print*,' cost_Geoid: egmunit              = ',egmunit
      print*,'             geoid_covariancefile = ',
     &                     geoid_covariancefile(ilo:ihi)
      print*,'             opening file ...'

      open(egmunit, file   = geoid_covariancefile(ilo:ihi),
     &              form   = 'unformatted',
     &              access = 'direct',
     &              recl   = ncshc*egmsize,
     &              status = 'old' )

c--   round up the quotient to get jsize
      jsize = ncshc / numberOfProcs
      if (jsize*numberOfProcs .lt. ncshc ) then
         jsize = jsize + 1
      end if

c--   compute sub-intervall of 1..ncshc for this processor
      joff = myProcId*jsize
      jbeg = 1 + joff
      jend = min( ncshc, joff+jsize )

c--   read part of geoid error covariance matrix
c--   and compute corresponding part of cost contribution
CADJ loop = parallel
      do j = jbeg,jend

         read(egmunit,rec=j) pinv

         factor = 0.
         do i = 1,ncShc
            factor = factor + pinv(i)*shc(i)*recip_rr
         enddo
         sphcost = sphcost + factor*shc(j)
      enddo

c--   close geoid error covariance file
      close(egmunit)

      _END_MASTER( mythid )

      call global_sum_r8( cost, mythid )
      sphcost = cost

#else  /* ALLOW_EGM96_ERROR_COV */
      sphcost = 0.

#endif /* ALLOW_EGM96_ERROR_COV */

      end


1	C $Header: /u/gcmpack/MITgcm/pkg/cost/Attic/cost_geoid.F,v 1.1.2.2 2002/04/04 10:58:58 heimbach Exp $
2
3	#include "COST_CPPOPTIONS.h"
4
5	subroutine cost_geoid(
6	O sphcost,
7	I shc,
8	I mythid
9	& )
10
11	c ==================================================================
12	c SUBROUTINE cost_geoid
13	c ==================================================================
14	c
15	c o Evaluate the cost function of the geoid contribution.
16	c
17	c started: Christian Eckert eckert@mit.edu 30-Jun-1999
18	c
19	c changed: Christian Eckert eckert@mit.edu 25-Feb-2000
20	c
21	c - Restructured the code in order to create a package
22	c for the MITgcmUV.
23	c
24	c changed: Ralf Giering Ralf.Giering@FastOpt.de 12-Jun-2001
25	c
26	c - totally rewrite for parallel processing
27	c
28	c ==================================================================
29	c SUBROUTINE cost_geoid
30	c ==================================================================
31
32	implicit none
33
34	c == global variables ==
35
36	#include "EEPARAMS.h"
37	#include "SIZE.h"
38	#include "PARAMS.h"
39
40	#include "ecco_cost.h"
41	#include "cost_sph.h"
42
43	c == routine arguments ==
44
45	_RL sphcost
46	_RL shc( ncShc )
47	integer mythid
48
49	#ifdef ALLOW_EGM96_ERROR_COV
50
51	c == local variables ==
52
53	Real*8 pinv( ncshc )
54
55	integer i,j
56	integer ilo,ihi
57	integer ilnblnk
58	integer ifnblnk
59	integer ireads
60	integer egmunit
61	integer egmsize
62	parameter( egmsize = 8 )
63
64	integer jsize, joff, jbeg, jend, j
65
66	_RL factor
67	_RL recip_rr
68
69	c == external functions ==
70	integer ilnblnk
71	external ilnblnk
72	integer ifnblnk
73	external ifnblnk
74
75	c == end of interface ==
76
77	c-- Only the master thread is doing I/O
78	_BEGIN_MASTER( mythid )
79
80	c-- initialise cost variable for all threads
81	c-- to allow usage of global_sum_r8
82	sphcost = 0. _d 0
83	do i=1,ncShc
84	pinv(i) = 0. _d 0
85	enddo
86
87	c-- Initialise variables.
88	recip_rr = recip_rsphere*recip_rsphere
89
90	c-- Open the geoid error covariance data file.
91	call mdsfindunit( egmunit, mythid )
92	ilo = ifnblnk(geoid_covariancefile)
93	ihi = ilnblnk(geoid_covariancefile)
94
95	print*
96	print*,' cost_Geoid: egmunit = ',egmunit
97	print*,' geoid_covariancefile = ',
98	& geoid_covariancefile(ilo:ihi)
99	print*,' opening file ...'
100
101	open(egmunit, file = geoid_covariancefile(ilo:ihi),
102	& form = 'unformatted',
103	& access = 'direct',
104	& recl = ncshc*egmsize,
105	& status = 'old' )
106
107	c-- round up the quotient to get jsize
108	jsize = ncshc / numberOfProcs
109	if (jsize*numberOfProcs .lt. ncshc ) then
110	jsize = jsize + 1
111	end if
112
113	c-- compute sub-intervall of 1..ncshc for this processor
114	joff = myProcId*jsize
115	jbeg = 1 + joff
116	jend = min( ncshc, joff+jsize )
117
118	c-- read part of geoid error covariance matrix
119	c-- and compute corresponding part of cost contribution
120	CADJ loop = parallel
121	do j = jbeg,jend
122
123	read(egmunit,rec=j) pinv
124
125	factor = 0.
126	do i = 1,ncShc
127	factor = factor + pinv(i)shc(i)recip_rr
128	enddo
129	sphcost = sphcost + factor*shc(j)
130	enddo
131
132	c-- close geoid error covariance file
133	close(egmunit)
134
135	_END_MASTER( mythid )
136
137	call global_sum_r8( cost, mythid )
138	sphcost = cost
139
140	#else /* ALLOW_EGM96_ERROR_COV */
141	sphcost = 0.
142
143	#endif /* ALLOW_EGM96_ERROR_COV */
144
145	end
146
147