Oasis3 4.0.2
|
00001 00002 !#define __VERBOSE 00003 #define __READ_BY_LOOKUP 00004 !---------------------------------------------------------------------- 00005 ! BOP 00006 ! 00007 ! !MODULE: mod_psmile_io 00008 ! !REMARKS: Initially programed by Reiner Vogelsang, SGI GmbH, and integrated 00009 ! into OASIS 3.0 and tested by Arnaud Caubel, CERFACS 00010 ! Bugs should be reported to reiner@sgi.com or Sophie.Valcke@cerfacs.fr. 00011 ! !REVISION HISTORY: 00012 ! 2003.01.31 Reiner Vogelsang Finalized initial version 00013 ! 2003.02.28 Arnaud Caubel Tests and correction of the READ_BY_LOOKUP 00014 ! 2003.04.24 Reiner Vogelsang Corrections with regard to CF naming convention. 00015 ! 2003.04.28 Reiner Vogelsang Added ProTex documentation header and 00016 ! static CVS versions string. 00017 ! 2004.03.08 Sophie Valcke Added the writting of bounds,lons and lats. 00018 ! 2004.04.07 Reiner Vogelsang Correction: bounds are written with 00019 ! number of corners as leading dimension 00020 ! 2004.19.06 Reiner Vogelsang Considering cyclic boundaries for longitudes. 00021 ! which lead to differences between corners 00022 ! and bounds by 360 degrees 00023 ! 2004.26.08 Reiner Vogelsang Correction of allocation scheme for arrays 00024 ! hosting the longitudes. 00025 ! 00026 ! !PUBLIC MEMBER FUNCTIONS: 00027 ! 00028 ! subroutine psmile_io_init_comp(id_error) 00029 ! This subroutine finds a local communicator for the MPI tasks 00030 ! of a model component participating in coupling. 00031 ! 00032 ! subroutine psmile_def_domains(id_error) 00033 ! This subroutine initializes the MPP_IO domains. 00034 ! 00035 ! subroutine psmile_def_files(id_error) 00036 ! This subroutine open files acccording to the kewords 00037 ! 00038 ! subroutine psmile_def_metadata(id_error) 00039 ! This subroutine writes CF metadata or finds a variable according 00040 ! to CF metadata information. 00041 ! 00042 ! subroutine psmile_close_files(id_error) 00043 ! This subroutine closes all opened psmile_io files 00044 ! 00045 ! subroutine psmile_io_cleanup(id_error) 00046 ! This subroutine deallocates all PSMILe I/O related data 00047 ! and exits MPP_IO 00048 ! 00049 ! subroutine psmile_read_8(id_port_id,rd_field,id_newtime) 00050 ! This subroutine reads REAL(KIND=4) data. 00051 ! 00052 ! subroutine psmile_read_4(id_port_id,rd_field,id_newtime) 00053 ! This subroutine reads REAL(KIND=8) data. 00054 ! 00055 ! subroutine psmile_write_4(id_port_id,rd_field,id_newtime) 00056 ! This subroutine writes REAL(KIND=4) data. 00057 ! 00058 ! subroutine psmile_write_8(id_port_id,rd_field,id_newtime) 00059 ! This subroutine writes REAL(KIND=8) data. 00060 ! 00061 ! subroutine indexi(n,arr,indx) 00062 ! This subroutine sorts by indexing. 00063 ! 00064 ! subroutine combine_with_date(cd_in,cd_mode,id_initial_date,cd_on) 00065 ! This subroutines appends the date in ISO format to a strings. 00066 ! 00067 00068 module mod_psmile_io 00069 #if !defined key_noIO 00070 ! !USES: 00071 use mod_kinds_model 00072 use mpp_mod_oa 00073 use mpp_io_mod_oa 00074 use mpp_domains_mod_oa 00075 ! !PUBLIC TYPES: 00076 !Highest unit number permittable by the OS. 00077 integer(kind=ip_intwp_p),parameter::PSMILE_MAX_UNIT=299 00078 !Highest unit number reserved by the OS. 00079 integer(kind=ip_intwp_p),parameter::PSMILE_MAX_RESERVED_UNIT=104 00080 !Lowest unit number reserved by the OS. 00081 integer(kind=ip_intwp_p),parameter::PSMILE_MIN_RESERVED_UNIT=100 00082 integer(kind=ip_intwp_p),parameter::PSMILE_MAX_CF_TABLE_ENTRIES=1024 00083 integer(kind=ip_intwp_p),allocatable::nbr_corners(:) 00084 00085 integer(kind=ip_intwp_p)::ig_mpp_io_comm,ig_color,ig_mpp_io_comm_size,ig_mpp_io_comm_rank 00086 integer(kind=ip_intwp_p)::ig_stackmax=1310720,ig_stackmaxd=655360 00087 integer(kind=ip_intwp_p)::ig_layout(2) 00088 !Global arrays to specify IO domains 00089 type(domain2D),allocatable::domain(:) 00090 type(domain1D),allocatable::xdom(:),ydom(:) 00091 logical,allocatable::lg_is_odd_distribution(:) 00092 !Global arrarys for opening files 00093 integer(kind=ip_intwp_p),allocatable::ig_units(:) 00094 integer(kind=ip_intwp_p),allocatable::ig_action(:) 00095 integer(kind=ip_intwp_p),allocatable::ig_form(:) 00096 integer(kind=ip_intwp_p),allocatable::ig_threading(:) 00097 integer(kind=ip_intwp_p),allocatable::ig_fileset(:) 00098 character(len=128),allocatable::cg_my_input_file(:) 00099 character(len=128),allocatable::cg_output_file(:) 00100 00101 00102 !Global arrays for netcdf metadata informations 00103 integer(kind=ip_intwp_p),allocatable::ig_var_index(:) 00104 integer(kind=ip_intwp_p),allocatable::ig_ntimes(:) 00105 !RV: 22.01.2003 00106 TYPE time_stamps 00107 ! real(kind=ip_double_p),pointer::rl_times(:) 00108 !AC:31/01 00109 real(kind=ip_double_p),dimension(:),pointer :: rl_times 00110 END TYPE time_stamps 00111 ! !PUBLIC DATA MEMBERS: 00112 type(time_stamps),allocatable :: dg_times(:) 00113 00114 type(fieldtype),allocatable::ig_vars(:),field(:),latij(:),lonij(:) 00115 type(fieldtype),allocatable::crnlatij(:), crnlonij(:) 00116 type(axistype),allocatable::x_axis(:),y_axis(:),z_axis(:),t_axis(:) 00117 type(axistype),allocatable::c_axis(:) 00118 ! These arrays are supposed to be part of a different module! 00119 ! However, no information is available in Oasis 3.0. So, let them live here 00120 ! for a while. 00121 !RV 00122 character(len=64),allocatable::cports_lgname(:) 00123 character(len=64),allocatable::cports_units(:) 00124 character(len=64),allocatable::cfldlab(:),cfunitslab(:) 00125 !RV 00126 ! !DESCRIPTION: 00127 ! This module contains the global data objects which hold informations for 00128 ! performing I/O on coupled fields and are the glue between routines defined 00129 ! within this file, mod_psmile_io.F90. 00130 ! !EOP 00131 !------------------------------------------------------------------------------- 00132 !$Id: mod_psmile_io.F90 2826 2010-12-10 11:14:21Z valcke $ 00133 !------------------------------------------------------------------------------- 00134 #endif 00135 end module mod_psmile_io 00136 00137 #if !defined key_noIO 00138 00139 subroutine psmile_io_init_comp(id_error) 00140 !-------------------------------------------------------------------------- 00141 00142 !Routine must be called in prism_init_comp. 00143 !The assumption is that each process of a model calls prism_init_comp. 00144 use mod_prism_proto 00145 use mod_comprism_proto 00146 use mod_psmile_io 00147 00148 implicit none 00149 #include <mpif.h> 00150 00151 !Arguments: 00152 !Error return code 00153 integer(kind=ip_intwp_p),intent(out)::id_error 00154 00155 !Local Variables: 00156 !Don't touch the following character variable! 00157 character(len=80),save::cl_version_string= 00158 '$Id: mod_psmile_io.F90 2826 2010-12-10 11:14:21Z valcke $' 00159 ! 00160 !Communicator which covers the whole model. 00161 integer(kind=ip_intwp_p)::il_local_comm 00162 00163 integer(kind=ip_intwp_p)::il_i,il_key 00164 integer(kind=ip_intwp_p),allocatable::il_cpl_ranks(:) 00165 integer(kind=ip_intwp_p)::ig_local_comm_grp,il_ncplprocs,ig_local_comm_size,il_mynum 00166 00167 integer(kind=ip_intwp_p)::il_file_unit,il_max_entry_id,il_no_of_entries 00168 integer(kind=ip_intwp_p)::il_pos 00169 character(len=64) :: cl_cfname,cl_cfunit 00170 logical :: ll_exist 00171 00172 00173 #ifdef __VERBOSE 00174 write(nulprt,*) ' ' 00175 write(nulprt,*) '| | | Entering - - psmile_io_init_comp' 00176 call flush(nulprt) 00177 #endif 00178 id_error=CLIM_OK 00179 il_local_comm=ig_local_comm 00180 00181 !Create a communicator which covers all processes of a model invoLved with 00182 !coupling. 00183 il_ncplprocs=kbcplproc(ig_mynummod) 00184 allocate(il_cpl_ranks(1:il_ncplprocs)) 00185 00186 do il_i=0,il_ncplprocs-1 00187 il_cpl_ranks(il_i+1)=il_i 00188 enddo 00189 00190 il_key=1 00191 call MPI_Comm_rank(il_local_comm,il_mynum,id_error) 00192 00193 ! 00194 !Create a color by finding the cpl ranks in the model communicator. 00195 ! 00196 if(ANY(il_cpl_ranks.eq.il_mynum)) then 00197 ig_color=1000*ig_mynummod 00198 else 00199 ig_color=1 00200 endif 00201 #ifdef __VERBOSE 00202 write(nulprt,*)'psmile_io_init_comp: rank ',il_mynum,' color ',ig_color 00203 #endif 00204 00205 ! 00206 !Do the splitting. 00207 ! 00208 #ifdef __VERBOSE 00209 write(nulprt,*) 'psmile_io_init_comp: Communicator splitting' 00210 #endif 00211 call MPI_Comm_split(il_local_comm,ig_color,il_key,ig_mpp_io_comm,id_error) 00212 call MPI_Comm_rank(ig_mpp_io_comm,ig_mpp_io_comm_rank,id_error) 00213 call MPI_Comm_size(ig_mpp_io_comm,ig_mpp_io_comm_size,id_error) 00214 00215 ! 00216 !Initialize FMS mpp_init,mpp_io,mpp_domains 00217 ! 00218 if (ig_color .eq. 1000*ig_mynummod ) then 00219 #ifdef __VERBOSE 00220 write(nulprt,*) 'psmile_io_init_comp: MPP_INIT : Model',ig_mynummod & 00221 ,ig_mpp_io_comm_rank 00222 #endif 00223 call mpp_init(mpp_comm=ig_mpp_io_comm,logfile=trim(cnaprt)) 00224 call mpp_domains_set_stack_size(ig_stackmaxd) 00225 !rv 00226 !rv Oasis 3 uses for the tracefiles the I/O units >= 1000. Since 00227 !rv mpp_io stores informations regarding to file specs in an array of 00228 !rv structs that array has to be allocated apropiately. 00229 !rv maxresunit specifies that the top maxresunits are excluded from the list 00230 !rv of the files to be closed on mpp_io_exit. 00231 !rv 00232 call mpp_io_init(maxunit=nulprt,maxresunit=25) 00233 call mpp_set_stack_size(ig_stackmax) 00234 endif 00235 00236 !After this point I assume that I/O routines are called by processes 00237 !involved with coupling only. 00238 00239 !RV At this point I am reading the CF name table definitions. 00240 00241 inquire(file='cf_name_table.txt',exist=ll_exist) 00242 if(ll_exist) then 00243 #ifdef __VERBOSE 00244 write(nulprt,*) 'psmile_io_init_comp: Reading CF name table!' 00245 call flush(nulprt) 00246 #endif 00247 00248 il_file_unit=nulprt+1 00249 open(file='cf_name_table.txt' & 00250 ,unit=il_file_unit & 00251 ,form='formatted' & 00252 ,status='old') 00253 00254 read(unit=il_file_unit,fmt=*,iostat=id_error) 00255 read(unit=il_file_unit,fmt=*,iostat=id_error) & 00256 il_max_entry_id, il_no_of_entries 00257 00258 if(id_error.ne.0) then 00259 write(nulprt,*) & 00260 'psmile_io_init_comp:cf_name_table.txt:' & 00261 ,' Reading of first record failed!' 00262 call flush(nulprt) 00263 call MPI_Abort(mpi_comm,0,mpi_error) 00264 endif 00265 00266 if(il_max_entry_id.gt.0.and. & 00267 il_max_entry_id.le.PSMILE_MAX_CF_TABLE_ENTRIES) then 00268 00269 allocate(cfldlab(1:il_max_entry_id),STAT=id_error) 00270 if(id_error.ne.0) then 00271 write(nulprt,*) & 00272 'psmile_io_init_comp: Allocation of cfldlab failed!' 00273 call flush(nulprt) 00274 call MPI_Abort(mpi_comm,0,mpi_error) 00275 endif 00276 00277 allocate(cfunitslab(1:il_max_entry_id),STAT=id_error) 00278 if(id_error.ne.0) then 00279 write(nulprt,*) & 00280 'psmile_io_init_comp: Allocation of cfunitslab failed!' 00281 call flush(nulprt) 00282 call MPI_Abort(mpi_comm,0,mpi_error) 00283 endif 00284 00285 else 00286 00287 write(nulprt,*) & 00288 'psmile_io_init_comp:cf_name_table.txt:' & 00289 ,'Your max. numlab is out of range !' 00290 call flush(nulprt) 00291 call MPI_Abort(mpi_comm,0,mpi_error) 00292 00293 endif 00294 00295 00296 read(unit=il_file_unit,fmt=*,iostat=id_error) 00297 do il_i=1,il_no_of_entries 00298 00299 read(unit=il_file_unit,fmt=*,iostat=id_error) & 00300 il_pos,cl_cfname,cl_cfunit 00301 00302 if(id_error.eq.0) then 00303 if(il_pos.le.il_max_entry_id) then 00304 cfldlab(il_pos)=trim(cl_cfname) 00305 cfunitslab(il_pos)=trim(cl_cfunit) 00306 else 00307 write(nulprt,*) & 00308 'psmile_io_init_comp:cf_name_table.txt:' & 00309 ,'Record ',il_i,': numlab =',il_pos,' out of range!' 00310 call flush(nulprt) 00311 call MPI_Abort(mpi_comm,0,mpi_error) 00312 endif 00313 else 00314 write(nulprt,*) & 00315 'psmile_io_init_comp:cf_name_table.txt:' & 00316 ,'Reading record ',il_i,' failed!' 00317 call flush(nulprt) 00318 call MPI_Abort(mpi_comm,0,mpi_error) 00319 endif 00320 00321 enddo 00322 00323 endif 00324 00325 00326 00327 !RV 00328 deallocate(il_cpl_ranks) 00329 #ifdef __VERBOSE 00330 write(nulprt,*) '| | | Leaving - - psmile_io_init_comp' 00331 write(nulprt,*) ' ' 00332 call flush(nulprt) 00333 #endif 00334 00335 return 00336 end 00337 00338 subroutine psmile_def_domains(id_error) 00339 !-------------------------------------------------------------------------- 00340 ! 00341 ! Routine defines the IO domain partitioning 00342 ! Called at the end of prism_enddef. 00343 ! 00344 use mod_prism_proto 00345 use mod_comprism_proto 00346 use mod_psmile_io 00347 use mod_psmile_io_interfaces, only:indexi 00348 implicit none 00349 00350 #include <mpif.h> 00351 00352 integer(kind=ip_intwp_p),intent(out)::id_error 00353 00354 integer(kind=ip_intwp_p)::il_port,il_i,il_j,il_ldx 00355 integer(kind=ip_intwp_p),allocatable::il_gextents(:,:,:) 00356 integer(kind=ip_intwp_p),allocatable::il_offsets(:,:,:) 00357 integer(kind=ip_intwp_p),allocatable::il_extentsx(:) 00358 integer(kind=ip_intwp_p),allocatable::il_extentsy(:) 00359 integer(kind=ip_intwp_p),allocatable::il_offsetx(:) 00360 integer(kind=ip_intwp_p),allocatable::il_offsety(:) 00361 logical,allocatable::ll_mask(:,:) 00362 integer(kind=ip_intwp_p),allocatable::il_pelist(:) 00363 integer(kind=ip_intwp_p)::il_ngx,il_ngy 00364 integer(kind=ip_intwp_p)::il_recvcount=1,il_sum 00365 00366 !Assumption: Processes beloging to the same model have the same 00367 ! strategy per port. Each process involved in the 00368 ! coupling has the same number of ports 00369 00370 #ifdef __VERBOSE 00371 write(nulprt,*) ' ' 00372 write(nulprt,*)'| | | Entering - - psmile_def_domains' 00373 call flush(nulprt) 00374 #endif 00375 00376 id_error=CLIM_ok 00377 00378 !Global arrays 00379 allocate(domain(1:nports)) 00380 allocate(xdom(1:nports)) 00381 allocate(ydom(1:nports)) 00382 allocate(lg_is_odd_distribution(1:nports)) 00383 00384 !Local arrays 00385 !il_gextents(:,1,:) carries extents in x, il_gextents(:,2,:) extents in y 00386 allocate(il_gextents(0:ig_mpp_io_comm_size-1,1:2,1:nports)) 00387 allocate(il_offsets(1:ig_mpp_io_comm_size,1:2,1:nports)) 00388 allocate(il_pelist(0:ig_mpp_io_comm_size-1)) 00389 allocate(il_extentsx(0:ig_mpp_io_comm_size-1)) 00390 allocate(il_extentsy(0:ig_mpp_io_comm_size-1)) 00391 allocate(il_offsetx(0:ig_mpp_io_comm_size-1)) 00392 allocate(il_offsety(0:ig_mpp_io_comm_size-1)) 00393 allocate(ll_mask(0:ig_mpp_io_comm_size-1,0:ig_mpp_io_comm_size-1)) 00394 00395 if(ig_color/ig_mynummod.ne.1000) then 00396 WRITE(nulprt,*)'psmile_def_domains: Called by proc which is not' & 00397 ,' involved with coupling!' 00398 call flush(nulprt) 00399 call MPI_ABORT(mpi_comm,0,mpi_err) 00400 endif 00401 00402 do il_port=1,nports 00403 #ifdef __VERBOSE 00404 WRITE(nulprt,*)'ig_def_state ',ig_def_state(il_port) 00405 #endif 00406 IF (ig_def_state(il_port) .eq. ip_ignout .or. & 00407 ig_def_state(il_port) .eq. ip_expout .or. & 00408 ig_def_state(il_port) .eq. ip_output .or. & 00409 ig_def_state(il_port) .eq. ip_input ) THEN 00410 #ifdef __VERBOSE 00411 WRITE(nulprt,*)'psmile_def_domains: port: ',il_port 00412 call flush(nulprt) 00413 #endif 00414 00415 if( mydist(CLIM_Strategy,il_port) .eq. CLIM_Serial) then 00416 00417 #ifdef __VERBOSE 00418 WRITE(nulprt,*)'psmile_def_domains: port: ',il_port,' CLIM_Serial' 00419 #endif 00420 00421 if(ig_mpp_io_comm_size.eq.1)then 00422 #ifdef __VERBOSE 00423 WRITE(nulprt,*)'psmile_def_domains: port: ',il_port & 00424 ,'(/1,mydist( CLIM_Segments+2,il_port),1,1/)' & 00425 ,1,mydist(CLIM_Segments+2,il_port),1,1 00426 ! call flush(nulprt) 00427 #endif 00428 ig_layout(1:2)=1 00429 call mpp_define_domains((/1,mydist( CLIM_Segments+2,il_port),1,1/)& 00430 , ig_layout,domain(il_port)) 00431 call mpp_get_domain_components( domain(il_port), xdom(il_port),& 00432 ydom(il_port) ) 00433 else 00434 write(nulprt, * ) 'psmile_def_domains: no. of procs for I/O ',& 00435 '> 1 and CLIM_Serial does not make sense! ' 00436 call MPI_ABORT(mpi_comm,0,mpi_err) 00437 00438 endif 00439 00440 else if( mydist(CLIM_Strategy,il_port) .eq. CLIM_Apple) then 00441 00442 #ifdef __VERBOSE 00443 !!$ WRITE(nulprt,*)'psmile_def_domains: port: ',il_port,' CLIM_Apple',mydist(:,il_port) 00444 #endif 00445 00446 !Assuming 1D decomposition. At the this point this is the only information 00447 !which I have considering CLIM_Apple! 00448 00449 ig_layout(1)=ig_mpp_io_comm_size 00450 ig_layout(2)=1 00451 00452 !Build a list of the local extents which exists on the procs involved with 00453 !coupling 00454 call MPI_Allgather(mydist(CLIM_Length+1,il_port),1,MPI_INTEGER,& 00455 il_gextents(0,1,il_port),il_recvcount,MPI_INTEGER,& 00456 ig_mpp_io_comm,id_error) 00457 call MPI_Allgather(mydist(CLIM_Offset+1,il_port),1,MPI_INTEGER,& 00458 il_offsets(1,1,il_port),il_recvcount,MPI_INTEGER,& 00459 ig_mpp_io_comm,id_error) 00460 00461 !Sort the extents and the PE-list according to the offsets 00462 00463 call indexi(ig_mpp_io_comm_size & 00464 ,il_offsets(1:ig_mpp_io_comm_size,1,il_port),il_pelist) 00465 00466 il_pelist(:)=il_pelist(:)-1 00467 il_extentsx(0:ig_mpp_io_comm_size-1)=& 00468 il_gextents(il_pelist(:),1,il_port) 00469 00470 !The max global index for the current 1D decomposition. 00471 il_sum=sum(il_extentsx(0:ig_mpp_io_comm_size-1)) 00472 00473 #ifdef __VERBOSE 00474 WRITE(nulprt,*)'psmile_def_domains: port: ',il_port & 00475 ,'(/1,il_sum,1,1/),pelist',1,il_sum,1,1,il_pelist & 00476 ,'extents',il_extentsx 00477 #endif 00478 !Define a pseudo 2D domain 00479 il_extentsy(0)=1 00480 call mpp_define_domains((/1,il_sum,1,1/),(/ig_mpp_io_comm_size,1/),& 00481 domain(il_port),& 00482 pelist=il_pelist,& 00483 xextent=il_extentsx,& 00484 yextent=il_extentsy(0:0)) 00485 call mpp_get_domain_components( domain(il_port), xdom(il_port),& 00486 ydom(il_port) ) 00487 00488 else if( mydist(CLIM_Strategy,il_port) .eq. CLIM_Box) then 00489 00490 #ifdef __VERBOSE 00491 WRITE(nulprt,*)'psmile_def_domains: port: ',il_port,' CLIM_Box' 00492 ! WRITE(nulprt,*)'psmile_def_domains: port: ',il_port,' CLIM_Box', & 00493 ! mydist(:,il_port),'::',myport(:,il_port) 00494 #endif 00495 !Build a list of the local extents which exists on the procs involved with 00496 !coupling 00497 !Block sizes in x 00498 call MPI_Allgather(mydist(CLIM_SizeX+1,il_port),1,MPI_INTEGER,& 00499 il_gextents(0,1,il_port),il_recvcount,MPI_INTEGER,& 00500 ig_mpp_io_comm,id_error) 00501 !Block sizes in y 00502 call MPI_Allgather(mydist(CLIM_Segments,il_port),1,MPI_INTEGER,& 00503 il_gextents(0,2,il_port),il_recvcount,MPI_INTEGER,& 00504 ig_mpp_io_comm,id_error) 00505 !Offsets of the blocks 00506 call MPI_Allgather(mydist(CLIM_Offset+1,il_port),1,MPI_INTEGER,& 00507 il_offsets(1,1,il_port),il_recvcount,MPI_INTEGER,& 00508 ig_mpp_io_comm,id_error) 00509 00510 call indexi(ig_mpp_io_comm_size & 00511 ,il_offsets(1:ig_mpp_io_comm_size,1,il_port),& 00512 il_pelist) 00513 il_pelist(:)=il_pelist(:)-1 00514 00515 il_extentsx(0:ig_mpp_io_comm_size-1)=& 00516 il_gextents(il_pelist(:),1,il_port) 00517 il_extentsy(0:ig_mpp_io_comm_size-1)=& 00518 il_gextents(il_pelist(:),2,il_port) 00519 #ifdef __VERBOSE 00520 ! write(nulprt,*)'psmile_def_domains:il_extentsx',il_extentsx 00521 ! write(nulprt,*)'psmile_def_domains:il_extentsy',il_extentsy 00522 #endif 00523 ! At this point I have to think about the domain decomposition 00524 ! 1st assumption : The leading dimension is a global parameter for 00525 ! the model. So, I don't do further checks on that. 00526 00527 il_ldx=mydist(CLIM_Segments+3,il_port)-mydist(CLIM_Offset+1,il_port) 00528 #ifdef __VERBOSE 00529 ! write(nulprt,*)'psmile_def_domains: LDX = ',il_ldx 00530 #endif 00531 00532 00533 !Count the number of offsets which are less than LDX to get the partitioning 00534 !in x-direction 00535 ig_layout(1)=& 00536 count(mask=il_offsets(:,1,il_port).lt.il_ldx) 00537 00538 if(ig_layout(1).gt.0) then 00539 ig_layout(2)=ig_mpp_io_comm_size/ig_layout(1) 00540 if(ig_layout(1)*ig_layout(2).ne.ig_mpp_io_comm_size) then 00541 lg_is_odd_distribution(il_port)=.true. 00542 else 00543 !check within a column y for different extents of y 00544 lg_is_odd_distribution(il_port)=.false. 00545 do il_j=0,ig_layout(2)-1 00546 do il_i=1,ig_layout(1)-1 00547 lg_is_odd_distribution(il_port)= & 00548 il_extentsy(ig_layout(1)*il_j).ne.& 00549 il_extentsy(il_i+ig_layout(1)*il_j) 00550 #ifdef __VERBOSE 00551 ! write(nulprt,*)'i j extentsy ',il_i,il_j, & 00552 ! il_extentsy(ig_layout(1)*il_j), & 00553 ! il_extentsy(il_i+ig_layout(1)*il_j) 00554 #endif 00555 enddo 00556 enddo 00557 00558 !check within a column x for different extents of x 00559 do il_i=0,ig_layout(1)-1 00560 do il_j=1,ig_layout(2)-1 00561 lg_is_odd_distribution(il_port)= & 00562 il_extentsx(il_i).ne.il_extentsx(il_i+ig_layout(1)*il_j) 00563 #ifdef __VERBOSE 00564 ! write(nulprt,*)'i j extentsy ',il_i,il_j, & 00565 ! il_extentsx(il_i), & 00566 ! il_extentsx(il_i+ig_layout(1)*il_j) 00567 #endif 00568 enddo 00569 enddo 00570 00571 endif 00572 00573 !Okay, fine! There are global divisors of the x- and y-axis! 00574 if(.not.lg_is_odd_distribution(il_port)) then 00575 00576 00577 il_ngx= il_ldx 00578 il_ngy=il_offsets(il_pelist(ig_mpp_io_comm_size-1)+1,1,il_port)& 00579 /il_ngx + il_extentsy(ig_mpp_io_comm_size-1) 00580 00581 #ifdef __VERBOSE 00582 WRITE(nulprt,*)'psmile_def_domains: port: ',il_port & 00583 ,'(/1,il_ngx,1,il_ngy/),pelist' & 00584 ,1,il_ngx,1,il_ngy,il_pelist 00585 #endif 00586 00587 call mpp_define_domains((/1,il_ngx,1,il_ngy/) & 00588 ,ig_layout,domain(il_port) & 00589 ,xextent=il_extentsx(0:ig_mpp_io_comm_size-1:ig_layout(2)) & 00590 ,yextent=il_extentsy(0:ig_mpp_io_comm_size-1:ig_layout(1)) & 00591 ,pelist=il_pelist) 00592 call mpp_get_domain_components( domain(il_port), xdom(il_port),& 00593 ydom(il_port) ) 00594 00595 else 00596 00597 !Okay, at this point it seems that have we have a block-structured 00598 !distribution. Need to calculate the offsets of the blocks in x- and y- 00599 !direction. These offsets migth be provided by prism_set_offset. 00600 ll_mask(:,:)=.false. 00601 do il_i=0,ig_mpp_io_comm_size-1 00602 ll_mask(il_pelist(il_i),il_pelist(il_i))=.true. 00603 enddo 00604 il_ngx=il_ldx 00605 il_offsetx(:)=il_offsets(il_pelist(0:ig_mpp_io_comm_size-1)+1,1,il_port) & 00606 -il_ngx* & 00607 (il_offsets(il_pelist(0:ig_mpp_io_comm_size-1)+1,1,il_port) & 00608 /il_ngx) 00609 il_offsety(:)=il_offsets(il_pelist(0:ig_mpp_io_comm_size-1)+1,1,il_port) & 00610 /il_ngx 00611 il_ngy=il_offsety(ig_mpp_io_comm_size-1) & 00612 +il_extentsy(ig_mpp_io_comm_size-1) 00613 ig_layout(1)=ig_mpp_io_comm_size 00614 ig_layout(2)=ig_mpp_io_comm_size 00615 00616 #ifdef __VERBOSE 00617 WRITE(nulprt,*)'psmile_def_domains: port: ',il_port & 00618 ,' :fall back: (/1,il_ngx,1,il_ngy/),pelist' & 00619 ,1,il_ngx,1,il_ngy,il_pelist,il_offsetx(:) & 00620 ,il_offsety(:) 00621 #endif 00622 00623 call mpp_define_domains((/1,il_ngx,1,il_ngy/) & 00624 ,ig_layout,domain(il_port) & 00625 ,xextent=il_extentsx & 00626 ,yextent=il_extentsy & 00627 ,pelist=il_pelist,maskmap=ll_mask & 00628 ,offsetx=il_offsetx,offsety=il_offsety) 00629 call mpp_get_domain_components( domain(il_port), xdom(il_port),& 00630 ydom(il_port) ) 00631 endif 00632 00633 else 00634 00635 write(nulprt, * ) 'psmile_def_domains: ig_layout(1) =0! ' 00636 call MPI_ABORT(mpi_comm,0,mpi_err) 00637 endif 00638 00639 else if( mydist(CLIM_Strategy,il_port) .eq. CLIM_Orange) then 00640 00641 write(nulprt, * ) 'psmile_def_domains: Strategy CLIM_orange',& 00642 'not yet implemented for I/O!' 00643 call MPI_ABORT(mpi_comm,0,mpi_err) 00644 00645 endif 00646 endif 00647 enddo 00648 00649 !Cleanup 00650 deallocate(il_pelist) 00651 deallocate(il_extentsx) 00652 deallocate(il_extentsy) 00653 deallocate(il_gextents) 00654 deallocate(il_offsets) 00655 deallocate(il_offsetx) 00656 deallocate(il_offsety) 00657 deallocate(ll_mask) 00658 00659 #ifdef __VERBOSE 00660 WRITE(nulprt,*)'| | | Leaving - - psmile_def_domains:' 00661 write(nulprt,*) ' ' 00662 call flush(nulprt) 00663 #endif 00664 00665 00666 end subroutine psmile_def_domains 00667 00668 subroutine psmile_def_files(id_error) 00669 !-------------------------------------------------------------------------- 00670 ! 00671 ! This subroutine open files acccording to the kewords 00672 ! OUTPUT: write by model issueing prism_put 00673 ! INPUT: read by issueing prism_get 00674 ! EXPOUT: write by model issueing prism_put and prims_get after receive. 00675 ! Usage for debugging data exchange model->Oasis->model 00676 ! IGNOUT: write by model issueing prism_put and prims_get after receive. 00677 ! Usage for debugging data exchange model->model 00678 ! 00679 ! Called at the end of prism_enddef. 00680 ! 00681 use mod_prism_proto 00682 use mod_comprism_proto 00683 use mod_psmile_io 00684 use mod_psmile_io_interfaces,only:combine_with_date 00685 use mod_psmile_date_and_time 00686 00687 implicit none 00688 00689 #include <mpif.h> 00690 00691 integer(kind=ip_intwp_p),intent(out)::id_error 00692 !Local variables 00693 integer(kind=ip_intwp_p)::il_i,il_port,il_unit 00694 !Test integer(kind=ip_intwp_p),allocatable::il_units(:) 00695 logical:: ll_opened,ll_exist 00696 00697 #ifdef __VERBOSE 00698 write(nulprt,*) ' ' 00699 write(nulprt,*)'| | | Entering - - psmile_def_files' 00700 call flush(nulprt) 00701 #endif 00702 00703 if(ig_color/ig_mynummod.ne.1000) then 00704 WRITE(nulprt,*)'psmile_def_files: Called by proc which is not' & 00705 ,' involved with coupling!' 00706 call flush(nulprt) 00707 call MPI_ABORT(mpi_comm,0,mpi_err) 00708 endif 00709 00710 if((.not.allocated(domain))) then 00711 WRITE(nulprt,*)'psmile_def_files:' & 00712 ,' Not called after psmile_def_domains!' 00713 call flush(nulprt) 00714 call MPI_ABORT(mpi_comm,0,mpi_err) 00715 endif 00716 00717 00718 id_error=CLIM_Ok 00719 00720 !Find unused units. 00721 !Problem is here that mpi codes usually are scattering/collecting from/on one PE 00722 ! which read/write. The rank of that PE is not known to PSMILE. 00723 !Algorithm: Let each PE involved with coupling do a search beginning from 00724 !PSMILE_MAX_UNIT. Check if the result is the same on each couple task. 00725 !Currently I don't know if the check is absolutely necessary. 00726 !Moreover, I assume that MPI is doing private I/O for each task, means 00727 !each opened file is handled individually per task by the OS . 00728 !I made my experience with UNICOS where file handling where global per default 00729 !running MPI in autotasking mode. In that case the user had to specify a 00730 !certain FFIO flag to get the normal behaviour. 00731 00732 !Test allocate(il_units(1:ig_mpp_io_comm_size)) 00733 allocate(ig_units(1:nports)) 00734 !RV:22.01.2003 00735 ig_units=-1 00736 00737 il_unit=PSMILE_MAX_UNIT 00738 do il_port=1,nports 00739 IF (ig_def_state(il_port) .eq. ip_ignout .or. & 00740 ig_def_state(il_port) .eq. ip_expout .or. & 00741 ig_def_state(il_port) .eq. ip_output .or. & 00742 ig_def_state(il_port) .eq. ip_input ) THEN 00743 do 00744 inquire(unit=il_unit,exist=ll_exist,opened=ll_opened) 00745 if( ll_exist .and. (.not.ll_opened)) exit 00746 il_unit=il_unit-1 00747 if(il_unit.eq.PSMILE_MAX_RESERVED_UNIT) & 00748 il_unit=PSMILE_MIN_RESERVED_UNIT-1 00749 if(il_unit.eq.7) then 00750 write(nulprt,*)'psmile_def_files: Could not find unit!' 00751 call MPI_ABORT(mpi_comm,0,mpi_err) 00752 endif 00753 enddo 00754 ig_units(il_port)=il_unit 00755 00756 #ifdef __VERBOSE 00757 write(nulprt,*)'psmile_def_files: Found FORTRAN I/O unit ', il_unit 00758 call flush(nulprt) 00759 #endif 00760 00761 il_unit=il_unit-1 00762 endif 00763 enddo 00764 00765 !Opening files! 00766 !I am opening the files such that read / write is performed on one file. 00767 !One can easily change that to distributed i/o changing the fileset attribute. 00768 00769 allocate(ig_action(1:nports)) 00770 allocate(ig_form(1:nports)) 00771 allocate(ig_threading(1:nports)) 00772 allocate(ig_fileset(1:nports)) 00773 allocate(cg_my_input_file(1:nports)) 00774 allocate(cg_output_file(1:nports)) 00775 00776 ! call_mpp_sync() 00777 do il_port=1,nports 00778 if(ig_def_state(il_port).eq.ip_input) then 00779 00780 ig_action(il_port)=MPP_RDONLY 00781 ig_form(il_port)=MPP_NETCDF 00782 ig_threading(il_port)=MPP_MULTI 00783 ig_fileset(il_port)=MPP_SINGLE 00784 00785 !I create here the name of the input file according to the convention 00786 !cports(il_port)_in.<year>-<month>-<day>T<hh>:<mm>:ss.nc 00787 !Arnaud, if you like you prefer approach of reading the name from the namcouple 00788 !initialize cg_my_input_file accordingly and comment out the call of 00789 !combine_with_date. 00790 00791 ! call combine_with_date(cports(il_port),'in',il_my_initial_date & 00792 ! ,cg_my_input_file(il_port)) 00793 00794 call mpp_open( ig_units(il_port),trim(cg_def_inpfile(il_port)) & 00795 , action=ig_action(il_port) & 00796 , form=ig_form(il_port) & 00797 , threading=ig_threading(il_port) & 00798 , fileset=ig_fileset(il_port)) 00799 00800 elseif(ig_def_state(il_port).eq.ip_output) then 00801 00802 ig_action(il_port)=MPP_OVERWR 00803 ig_form(il_port)=MPP_NETCDF 00804 ig_threading(il_port)=MPP_SINGLE 00805 00806 call combine_with_date(cports(il_port),'out',ig_inidate & 00807 ,cg_output_file(il_port)) 00808 00809 call mpp_open( ig_units(il_port), trim(cg_output_file(il_port)) & 00810 , action=ig_action(il_port) & 00811 , form=ig_form(il_port) & 00812 , threading=ig_threading(il_port)) 00813 00814 00815 elseif(ig_def_state(il_port).eq.ip_expout) then 00816 00817 ig_action(il_port)=MPP_OVERWR 00818 ig_form(il_port)=MPP_NETCDF 00819 ig_threading(il_port)=MPP_SINGLE 00820 00821 call combine_with_date(cports(il_port),'out',ig_inidate & 00822 ,cg_output_file(il_port)) 00823 00824 call mpp_open( ig_units(il_port),trim(cg_output_file(il_port)) & 00825 , action=ig_action(il_port) & 00826 , form=ig_form(il_port) & 00827 , threading=ig_threading(il_port)) 00828 00829 elseif(ig_def_state(il_port).eq.ip_ignout) then 00830 00831 ig_action(il_port)=MPP_OVERWR 00832 ig_form(il_port)=MPP_NETCDF 00833 ig_threading(il_port)=MPP_SINGLE 00834 00835 call combine_with_date(cg_ignout_field(il_port),'out',ig_inidate & 00836 ,cg_output_file(il_port)) 00837 00838 call mpp_open( ig_units(il_port), trim(cg_output_file(il_port)) & 00839 , action=ig_action(il_port) & 00840 , form=ig_form(il_port) & 00841 , threading=ig_threading(il_port)) 00842 endif 00843 enddo 00844 00845 00846 #ifdef __VERBOSE 00847 write(nulprt,*)'| | | Leaving - - psmile_def_files' 00848 write(nulprt,*) ' ' 00849 call flush(nulprt) 00850 #endif 00851 end subroutine psmile_def_files 00852 00853 subroutine psmile_def_metadata(id_error) 00854 !-------------------------------------------------------------------------- 00855 ! 00856 !This routine defines the metadata headers of netcdf files case or writing. 00857 !In case or reading it defines the index to lookup a certain variable in 00858 !a netcdf file 00859 ! Called at the end of prism_enddef. 00860 use mod_kinds_model 00861 use mod_psmile_io_interfaces,only:combine_with_date 00862 use mod_prism_proto 00863 use mod_comprism_proto 00864 use mod_psmile_io 00865 implicit none 00866 00867 #include <mpif.h> 00868 00869 integer(kind=ip_intwp_p),intent(out)::id_error 00870 !Local 00871 character(len=64)::cl_varname 00872 integer(kind=ip_intwp_p)::il_port,il_i 00873 integer(kind=ip_intwp_p)::il_xbegin,il_xend,il_ybegin,il_yend 00874 integer(kind=ip_intwp_p)::il_ndim,il_nvar,il_natt,il_ntime 00875 integer(kind=ip_intwp_p)::il_crn,il_ii,il_jj 00876 00877 type(fieldtype),allocatable::il_vars(:) 00878 ! 00879 ! Dimensions of the lons and lats contained in the OASIS 3 grid file 00880 integer(kind=ip_intwp_p)::lens_lat(4),lens_lon(4) 00881 INTEGER(kind=ip_intwp_p)::llo1_tmp, llo2_tmp 00882 INTEGER(kind=ip_intwp_p)::lla1_tmp, lla2_tmp 00883 REAL(kind=ip_double_p), ALLOCATABLE :: lon(:,:),lat(:,:) 00884 REAL(kind=ip_double_p), ALLOCATABLE :: loncrn(:,:,:),latcrn(:,:,:) 00885 LOGICAL :: lg_exist 00886 CHARACTER(len=8):: cg_clim_cgrdnamnc, clstrg 00887 INTEGER(kind=ip_intwp_p):: nc_grdid, il_varid, il_x, il_y 00888 REAL(kind=ip_double_p), ALLOCATABLE :: rl_lons(:,:), rl_lats(:,:) 00889 REAL(kind=ip_double_p), ALLOCATABLE :: rl_lonscrn(:,:,:), rl_latscrn(:,:,:) 00890 REAL(kind=ip_double_p):: difflon 00891 REAL(kind=ip_double_p):: threesixty 00892 LOGICAL :: llylat,llxlon, ll_notfound 00893 00894 #ifdef __VERBOSE 00895 write(nulprt,*)' ' 00896 write(nulprt,*)'| | | Entering -- psmile_def_metada' 00897 call flush(nulprt) 00898 #endif 00899 00900 if(ig_color/ig_mynummod.ne.1000) then 00901 WRITE(nulprt,*)'psmile_def_metdata: Called by proc which is not' & 00902 ,' involved with coupling!' 00903 call flush(nulprt) 00904 call MPI_ABORT(mpi_comm,0,mpi_err) 00905 endif 00906 00907 if((.not.allocated(ig_units)).and.(.not.allocated(domain))) then 00908 WRITE(nulprt,*)'psmile_def_metdata:' & 00909 ,' Not called after psmile_def_domains' & 00910 ,' and psmile_def_files!' 00911 call flush(nulprt) 00912 call MPI_ABORT(mpi_comm,0,mpi_err) 00913 endif 00914 00915 00916 id_error=CLIM_Ok 00917 threesixty=360.0 00918 00919 allocate(ig_var_index(nports)) 00920 allocate(ig_vars(1:nports)) 00921 ! 2005-07-20, Luis Kornblueh, MPIMet 00922 ! for later association test the pointer must be initialized, otherwise 00923 ! the program might crash - that is a very common problem! mpp_io is providing 00924 ! an initial axis which is extended by the nullifyd data pointer. 00925 ig_vars(:) = default_field 00926 allocate(ig_ntimes(1:nports)) 00927 00928 !RV:22.01.2003 00929 #ifdef __READ_BY_LOOKUP 00930 allocate(dg_times(1:nports)) 00931 !AC:31/01 00932 DO il_i = 1,nports 00933 nullify(dg_times(il_i)%rl_times) 00934 ENDDO 00935 #endif 00936 00937 allocate(x_axis(1:nports)) 00938 ! 2005-07-20, Luis Kornblueh, MPIMet 00939 x_axis(:) = default_axis 00940 allocate(y_axis(1:nports)) 00941 ! 2005-07-20, Luis Kornblueh, MPIMet 00942 y_axis(:) = default_axis 00943 allocate(c_axis(1:nports)) 00944 ! 2005-07-20, Luis Kornblueh, MPIMet 00945 c_axis(:) = default_axis 00946 #ifdef __VERBOSE 00947 write(nulprt,*)'psmile_def_metada: Number of corners: ',ig_noc 00948 call flush(nulprt) 00949 #endif 00950 allocate(nbr_corners(1:nports)) 00951 nbr_corners=ig_noc 00952 !Was not allocated. Oasis sends 2D arrays. In case of bundles or 3d arrays 00953 !we have to declare z_axis + plus one for bundles. 00954 !Left here as a reminder! 00955 ! allocate(z_axis(1:nports)) 00956 ! 2005-07-20, Luis Kornblueh, MPIMet 00957 ! z_axis(:) = default_axis 00958 allocate(t_axis(1:nports)) 00959 ! 2005-07-20, Luis Kornblueh, MPIMet 00960 t_axis(:) = default_axis 00961 allocate(field(1:nports)) 00962 ! 2005-07-20, Luis Kornblueh, MPIMet 00963 field(:) = default_field 00964 allocate(lonij(1:nports)) 00965 ! 2005-07-20, Luis Kornblueh, MPIMet 00966 lonij(:) = default_field 00967 allocate(latij(1:nports)) 00968 ! 2005-07-20, Luis Kornblueh, MPIMet 00969 latij(:) = default_field 00970 allocate(crnlonij(1:nports)) 00971 ! 2005-07-20, Luis Kornblueh, MPIMet 00972 crnlonij(:) = default_field 00973 allocate(crnlatij(1:nports)) 00974 ! 2005-07-20, Luis Kornblueh, MPIMet 00975 crnlatij(:) = default_field 00976 00977 !The following files are allocated here since no informations 00978 !are given in Oasis regarding to long names and units. I am filling them with 00979 !dummy informations. I put the data structure into mod_psmile_io. 00980 allocate(cports_lgname(1:nports)) 00981 allocate(cports_units(1:nports)) 00982 cports_lgname='BLANK' 00983 cports_units='1' 00984 00985 00986 DO il_port=1,nports 00987 IF(ig_def_state(il_port).EQ.ip_input) THEN 00988 00989 !Get information with regard to no. of dimensions, variables, attributes 00990 !and time stamps. 00991 call mpp_get_info(ig_units(il_port)& 00992 ,il_ndim,il_nvar,il_natt,il_ntime) 00993 00994 ig_ntimes(il_port)=il_ntime 00995 00996 allocate(il_vars(1:il_nvar)) 00997 !RV,16.12.2002:Get the all ids of the fields within the file. 00998 call mpp_get_fields (ig_units(il_port),il_vars(:)) 00999 01000 !Loop over variables and find the index of cports(il_port). 01001 ig_var_index(il_port)=-1 01002 do il_i=1,il_nvar 01003 01004 call mpp_get_atts(il_vars(il_i),name=cl_varname) 01005 if(trim(cl_varname).eq.trim(cports(il_port)) ) then 01006 ig_var_index(il_port)=il_i 01007 ig_vars(il_port)=il_vars(il_i) 01008 endif 01009 01010 enddo 01011 01012 !Check: Did I find the variable? If not, the user has given the wrong 01013 ! input file. 01014 if(ig_var_index(il_port).lt.0) then 01015 01016 write(nulprt, * ) & 01017 'psmile_def_metadata: Variable name given in namcouple '& 01018 ,cports(il_port) & 01019 ,'could not be found in file ',cg_my_input_file(il_port) & 01020 ,' port number: ',il_port 01021 call MPI_ABORT(mpi_comm,0,mpi_err) 01022 01023 endif 01024 01025 deallocate(il_vars) 01026 01027 elseif((ig_def_state(il_port).eq.ip_output).or. & 01028 (ig_def_state(il_port).eq.ip_expout).or. & 01029 (ig_def_state(il_port).eq.ip_ignout))then 01030 01031 !Define the global axis X and Y. 01032 !I have no informations on the units of the axis. I assume that 01033 !the axis are indices in order to perform a mapping like lat(i,j). 01034 call mpp_get_global_domain(domain(il_port) & 01035 ,xbegin=il_xbegin,xend=il_xend & 01036 ,ybegin=il_ybegin,yend=il_yend) 01037 call mpp_write_meta(ig_units(il_port) & 01038 ,'Conventions' & 01039 ,cval='CF-1.0') 01040 #ifdef __VERBOSE 01041 write(nulprt,*)'psmile_def_metadata: x-axis ',il_xbegin,il_xend 01042 #endif 01043 01044 call mpp_write_meta(ig_units(il_port),x_axis(il_port) & 01045 ,'X','1','global index space for x' & 01046 ,domain=xdom(il_port) & 01047 ,data=(/(il_i-1.,il_i=il_xbegin,il_xend)/) ) 01048 01049 #ifdef __VERBOSE 01050 write(nulprt,*)'psmile_def_metadata: y-axis ',il_ybegin,il_yend 01051 call flush(nulprt) 01052 #endif 01053 01054 call mpp_write_meta(ig_units(il_port),y_axis(il_port) & 01055 ,'Y','1','global index space for y' & 01056 ,domain=ydom(il_port) & 01057 ,data=(/(il_i-1.,il_i=il_ybegin,il_yend)/) ) 01058 #ifdef __VERBOSE 01059 write(nulprt,*)'psmile_def_metadata: c-axis ' 01060 #endif 01061 call mpp_write_meta(ig_units(il_port),c_axis(il_port) & 01062 ,'C','1','global index space for c' & 01063 ,data=(/(il_i-1.,il_i=1,nbr_corners(il_port))/) ) 01064 01065 01066 !Define the unlimited time axis 01067 call combine_with_date('second','since',ig_inidate & 01068 ,cl_varname) 01069 01070 call mpp_write_meta( ig_units(il_port), t_axis(il_port) & 01071 , 'time', trim(cl_varname), 'Time' ) 01072 !Define the real coordinates , here lattitudes and longitudes. 01073 ! I have no information with regard to latitude and longitude 01074 !RV:17.04.2003 Wright this information anyway to pass the CF check. 01075 call mpp_write_meta( ig_units(il_port), latij(il_port) & 01076 , (/x_axis(il_port),y_axis(il_port)/) & 01077 , 'lat', 'degrees_north' & 01078 , 'latitude') 01079 call mpp_write_meta( ig_units(il_port), lonij(il_port) & 01080 , (/x_axis(il_port),y_axis(il_port)/) & 01081 , 'lon', 'degrees_east' & 01082 , 'longitude') 01083 !rv,sgi Number of corners has to be the leading dimension according the 01084 !rv,sgi CF convention 01085 call mpp_write_meta( ig_units(il_port), crnlatij(il_port) & 01086 , (/c_axis(il_port),x_axis(il_port),y_axis(il_port)/) & 01087 , 'bounds_lat', 'degrees_north' & 01088 , 'corner latitudes') 01089 call mpp_write_meta( ig_units(il_port), crnlonij(il_port) & 01090 , (/c_axis(il_port),x_axis(il_port),y_axis(il_port)/) & 01091 , 'bounds_lon', 'degrees_east' & 01092 , 'corner longitudes') 01093 !Define that variable cports(il_port) of unit cports_units(il_port) and 01094 !longname cports_lgname(il_port) will be written as field field(il_port) 01095 !according to the three axis x_axis(il_port),y_axis(il_port) and 01096 !t_axis(il_port). The precision is double (pack=1). 01097 !RV 01098 ! Initialize cports_lgname and cports_units according to the CF name table 01099 ! exploiting the field label. 01100 ! 01101 IF(allocated(cfldlab)) then 01102 IF(ig_def_numlab(il_port).gt.0.and. & 01103 ig_def_numlab(il_port).le.PSMILE_MAX_CF_TABLE_ENTRIES) then 01104 cports_lgname(il_port)=cfldlab(ig_def_numlab(il_port)) 01105 cports_units(il_port)=cfunitslab(ig_def_numlab(il_port)) 01106 ENDIF 01107 ENDIF 01108 !RV 01109 IF (ig_def_state(il_port).EQ.ip_ignout) THEN 01110 call mpp_write_meta( ig_units(il_port), field(il_port) & 01111 , (/x_axis(il_port),y_axis(il_port) & 01112 ,t_axis(il_port)/) & 01113 ,trim(cg_ignout_field(il_port)) & 01114 ,trim(cports_units(il_port)) & 01115 , trim(cports_lgname(il_port)), pack=1 ) 01116 #ifdef __DEBUG 01117 write(nulprt,*)'psmile_def_metadata 1:',ig_def_numlab(il_port) & 01118 ,' ',trim(cg_ignout_field(il_port)) & 01119 ,' ',trim(cports_lgname(il_port)) & 01120 ,' ',trim(cports_units(il_port)) 01121 call flush(nulprt) 01122 #endif 01123 ELSE 01124 call mpp_write_meta( ig_units(il_port), field(il_port) & 01125 , (/x_axis(il_port),y_axis(il_port) & 01126 ,t_axis(il_port)/) & 01127 ,trim(cports(il_port)) & 01128 ,trim(cports_units(il_port)) & 01129 , trim(cports_lgname(il_port)), pack=1 ) 01130 #ifdef __DEBUG 01131 write(nulprt,*)'psmile_def_metadata 2:',ig_def_numlab(il_port) & 01132 ,' ', trim(cports(il_port)) & 01133 ,' ',trim(cports_lgname(il_port)) & 01134 ,' ',trim(cports_units(il_port)) 01135 call flush(nulprt) 01136 #endif 01137 ENDIF 01138 01139 ! 01140 ! Define an additional attribute bounds for lon and lats. 01141 ! 01142 call mpp_write_meta(ig_units(il_port) & 01143 ,mpp_get_id(latij(il_port)) & 01144 , 'bounds',cval='bounds_lat') 01145 call mpp_write_meta(ig_units(il_port) & 01146 ,mpp_get_id(lonij(il_port)) & 01147 , 'bounds',cval='bounds_lon') 01148 !Define an additional attribute that field(il_port) lives on the coordinates 01149 !lon and lat. This is a CF requirement! 01150 01151 call mpp_write_meta(ig_units(il_port) & 01152 ,mpp_get_id(field(il_port)) & 01153 , 'coordinates',cval='lon lat') 01154 01155 !Flush the metadata information to the file connected to unit 01156 !ig_units(il_port). 01157 01158 call mpp_write(ig_units(il_port), x_axis(il_port) ) 01159 call mpp_write(ig_units(il_port), y_axis(il_port) ) 01160 ! 01161 !rv Metadata should be flushed from the struct c_axis 01162 !rv into the NetCDF file 01163 ! 01164 call mpp_write(ig_units(il_port), c_axis(il_port) ) 01165 ! 01166 IF (mpp_pe() .EQ. mpp_root_pe()) THEN 01167 WRITE(nulprt,*)'mod_psmile_io - - check grids.nc file:' 01168 cg_clim_cgrdnamnc=cg_clim_cgrdnam//'.nc' 01169 INQUIRE(FILE = cg_clim_cgrdnamnc, EXIST = lg_exist) 01170 IF (.NOT. lg_exist) THEN 01171 WRITE(nulprt,*) ' - no grids.nc file found' 01172 WRITE(nulprt,*) 'lats and lons will not be written to output files' 01173 CALL flush(nulprt) 01174 ELSE 01175 WRITE(nulprt,*)'mod_psmile_io - - grids.nc file found' 01176 ! Read the lats and the lons 01177 CALL mpp_open(nc_grdid, 'grids.nc', & 01178 action=MPP_RDONLY, & 01179 form=MPP_NETCDF, & 01180 threading=MPP_SINGLE) 01181 CALL mpp_get_info(nc_grdid, il_ndim, il_nvar, & 01182 il_natt, il_ntime ) 01183 ALLOCATE(il_vars(1:il_nvar)) 01184 CALL mpp_get_fields (nc_grdid,il_vars(:)) 01185 ! 01186 ! Longitudes 01187 ! 01188 clstrg = cga_clim_locator(il_port)//cg_clim_lonsuf 01189 #ifdef __VERBOSE 01190 WRITE(nulprt,*) 'clstrg = ', clstrg 01191 #endif 01192 il_i = 0 01193 ll_notfound = .TRUE. 01194 DO WHILE (ll_notfound) 01195 il_i = il_i + 1 01196 IF (il_i .GT. il_nvar) THEN 01197 WRITE(nulprt,*)'psmile_def_metdata: Pb in file grids.nc: No longitudes' & 01198 ,'for grid', cga_clim_locator(il_port) 01199 WRITE(nulprt,*)'Calling MPI_ABORT ' 01200 CALL flush(nulprt) 01201 CALL MPI_ABORT(mpi_comm,0,mpi_err) 01202 ENDIF 01203 CALL mpp_get_atts(il_vars(il_i),name=cl_varname,siz=lens_lon) 01204 #ifdef __VERBOSE 01205 WRITE(nulprt,*) 'cl_varname = ', cl_varname 01206 #endif 01207 ! 01208 !RV,sgi ALLOCATE(lon(1:il_xend-il_xbegin+1,1:il_yend-il_ybegin+1)) 01209 ! 01210 IF(TRIM(cl_varname).EQ.TRIM(clstrg) ) THEN 01211 ll_notfound = .FALSE. 01212 #ifdef __VERBOSE 01213 write(nulprt,*)'psmile_def_metadata: cl_varname= ' & 01214 ,trim(cl_varname) 01215 call flush(nulprt) 01216 #endif 01217 IF(.NOT.ALLOCATED(lon)) & 01218 ALLOCATE(lon(1:il_xend-il_xbegin+1,1:il_yend-il_ybegin+1)) 01219 ALLOCATE(rl_lons(lens_lon(1),lens_lon(2))) 01220 01221 !NOTE: I am using tindex = 1 in order to circumvent a 01222 !peculiarity of mpp_io: grids.nc contains the lons/lats 01223 !declared as lon(i,j)/lat(i,j). mpp_io interprets i and j as 01224 !axes. However, the axis i and j have no associated data. 01225 !Therefore, i and j are NOT taken as spatial axes but 01226 !as something related to a time stamp. Until that "feature" 01227 !is corrected use the call below. 01228 ! 01229 CALL mpp_read(nc_grdid,il_vars(il_i),rl_lons,tindex=1) 01230 ! 01231 !rv,sgi< 01232 !If the user has specified INVERT/REVERSE and 01233 !$CORLAT=NORSUD and/or $CORLON=ESTWST in the 01234 !preprocessing section of the field we INVERT the 01235 !ordering of longitudes such that data and grid are 01236 !are matching. The assertion is that 01237 !grid informations are initialized SUDNOR and WSTEST. 01238 01239 if((ig_def_invert(il_port).gt.1).or. & 01240 (ig_def_reverse(il_port).gt.1)) then 01241 llxlon=ig_def_invert(il_port).eq.3 & 01242 .or.ig_def_invert(il_port).eq.4 & 01243 .or.ig_def_reverse(il_port).eq.3 & 01244 .or.ig_def_reverse(il_port).eq.4 01245 llylat=ig_def_invert(il_port).eq.2 & 01246 .or.ig_def_invert(il_port).eq.4 & 01247 .or.ig_def_reverse(il_port).eq.2 & 01248 .or.ig_def_reverse(il_port).eq.4 01249 call inv2d (rl_lons,lens_lon(1),lens_lon(2) & 01250 ,llxlon,llylat) 01251 endif 01252 !rv,sgi> 01253 !Reshape the read array such that it matches the current 01254 !axes of the output file. 01255 ! 01256 lon=RESHAPE(rl_lons,shape=SHAPE(lon)) 01257 ! 01258 !Dump it ot the output file. 01259 ! 01260 CALL mpp_write(ig_units(il_port),lonij(il_port),lon) 01261 DEALLOCATE(rl_lons) 01262 ENDIF 01263 ENDDO 01264 ! 01265 !Latitudes 01266 ! 01267 clstrg = cga_clim_locator(il_port)//cg_clim_latsuf 01268 #ifdef __VERBOSE 01269 WRITE(nulprt,*) 'clstrg = ', clstrg 01270 #endif 01271 il_i = 0 01272 ll_notfound = .TRUE. 01273 DO WHILE (ll_notfound) 01274 il_i = il_i + 1 01275 IF (il_i .GT. il_nvar) THEN 01276 WRITE(nulprt,*)'psmile_def_metdata: Pb in file grids.nc: No latitudes' & 01277 ,'for grid', cga_clim_locator(il_port) 01278 WRITE(nulprt,*)'Calling MPI_ABORT ' 01279 CALL flush(nulprt) 01280 CALL MPI_ABORT(mpi_comm,0,mpi_err) 01281 ENDIF 01282 ! 01283 CALL mpp_get_atts(il_vars(il_i),name=cl_varname,siz=lens_lat) 01284 #ifdef __VERBOSE 01285 WRITE(nulprt,*) 'cl_varname = ', cl_varname 01286 #endif 01287 IF(TRIM(cl_varname).EQ.TRIM(clstrg) ) THEN 01288 ll_notfound = .FALSE. 01289 #ifdef __VERBOSE 01290 WRITE(nulprt,*)'psmile_def_metadata: cl_varname= ' & 01291 ,trim(cl_varname) 01292 CALL flush(nulprt) 01293 #endif 01294 ALLOCATE(lat(1:il_xend-il_xbegin+1,1:il_yend-il_ybegin+1)) 01295 ALLOCATE(rl_lats(lens_lat(1),lens_lat(2))) 01296 01297 CALL mpp_read(nc_grdid,il_vars(il_i),rl_lats,tindex=1) 01298 !rv,sgi< 01299 if((ig_def_invert(il_port).gt.1).or. & 01300 (ig_def_reverse(il_port).gt.1)) then 01301 llxlon=ig_def_invert(il_port).eq.3 & 01302 .or.ig_def_invert(il_port).eq.4 & 01303 .or.ig_def_reverse(il_port).eq.3 & 01304 .or.ig_def_reverse(il_port).eq.4 01305 llylat=ig_def_invert(il_port).eq.2 & 01306 .or.ig_def_invert(il_port).eq.4 & 01307 .or.ig_def_reverse(il_port).eq.2 & 01308 .or.ig_def_reverse(il_port).eq.4 01309 call inv2d (rl_lats,lens_lat(1),lens_lat(2) & 01310 ,llxlon,llylat) 01311 endif 01312 !rv,sgi> 01313 lat=reshape(rl_lats,shape=shape(lat)) 01314 call mpp_write(ig_units(il_port),latij(il_port),lat) 01315 DEALLOCATE(rl_lats) 01316 ENDIF 01317 ENDDO 01318 ! 01319 ! Longitude corners 01320 ! 01321 clstrg = cga_clim_locator(il_port)//crn_clim_lonsuf 01322 #ifdef __VERBOSE 01323 WRITE(nulprt,*) 'clstrg = ', clstrg 01324 #endif 01325 ll_notfound = .TRUE. 01326 01327 DO il_i=1,il_nvar 01328 01329 CALL mpp_get_atts(il_vars(il_i),name=cl_varname,siz=lens_lon) 01330 #ifdef __VERBOSE 01331 WRITE(nulprt,*) 'cl_varname = ', cl_varname 01332 #endif 01333 IF(TRIM(cl_varname).EQ.TRIM(clstrg) ) THEN 01334 ll_notfound = .FALSE. 01335 WRITE(nulprt,*) 'cl_varname = clstrg Longitude corners' 01336 #ifdef __VERBOSE 01337 write(nulprt,*) & 01338 'psmile_def_metada, lon, clstrg: ',clstrg,lens_lon 01339 #endif 01340 nbr_corners(il_port)=lens_lon(3) 01341 llo1_tmp = lens_lon(1) 01342 llo2_tmp = lens_lon(2) 01343 ALLOCATE(rl_lonscrn(lens_lon(1),lens_lon(2) & 01344 ,nbr_corners(il_port))) 01345 ALLOCATE(loncrn(nbr_corners(il_port) & 01346 ,1:il_xend-il_xbegin+1 & 01347 ,1:il_yend-il_ybegin+1)) 01348 lens_lon(1)=il_xend-il_xbegin+1 01349 lens_lon(2)=il_yend-il_ybegin+1 01350 lens_lon(3)=nbr_corners(il_port) 01351 ! 01352 CALL mpp_read(nc_grdid,il_vars(il_i),rl_lonscrn,tindex=1) 01353 !rv,sgi< 01354 if((ig_def_invert(il_port).gt.1).or. & 01355 (ig_def_reverse(il_port).gt.1)) then 01356 01357 llxlon=ig_def_invert(il_port).eq.3 & 01358 .or.ig_def_invert(il_port).eq.4 & 01359 .or.ig_def_reverse(il_port).eq.3 & 01360 .or.ig_def_reverse(il_port).eq.4 01361 llylat=ig_def_invert(il_port).eq.2 & 01362 .or.ig_def_invert(il_port).eq.4 & 01363 .or.ig_def_reverse(il_port).eq.2 & 01364 .or.ig_def_reverse(il_port).eq.4 01365 call inv2dcrn (rl_lonscrn,llo1_tmp,llo2_tmp & 01366 ,nbr_corners(il_port),llxlon,llylat) 01367 endif 01368 !rv,sgi> 01369 01370 ! 01371 !Reshape the read array such that it matches the current 01372 !axes of the output file. 01373 ! 01374 call transp_crn2cf(rl_lonscrn,loncrn,lens_lon) 01375 !rv,sgi< 01376 ! 01377 ! Check for cyclic boundaries ! 01378 ! For plot routines using a mesh-fill method 01379 ! the distances of the cell centers to the corners 01380 ! have to be less equal 360. degrees. 01381 ! 01382 do il_crn=1,nbr_corners(il_port) 01383 do il_jj=1,il_yend-il_ybegin+1 01384 do il_ii=1,il_xend-il_xbegin+1 01385 difflon=lon(il_ii,il_jj)- & 01386 loncrn(il_crn,il_ii,il_jj) 01387 if(difflon.gt.300.0) then 01388 loncrn(il_crn,il_ii,il_jj)= & 01389 loncrn(il_crn,il_ii,il_jj)+threesixty 01390 elseif(difflon.lt.-300.0) then 01391 loncrn(il_crn,il_ii,il_jj)= & 01392 loncrn(il_crn,il_ii,il_jj)-threesixty 01393 endif 01394 enddo 01395 enddo 01396 enddo 01397 !rv,sgi> 01398 ! 01399 !Dump it to the output file. 01400 ! 01401 CALL mpp_write(ig_units(il_port),crnlonij(il_port),loncrn) 01402 DEALLOCATE(rl_lonscrn) 01403 DEALLOCATE(loncrn) 01404 ENDIF 01405 ENDDO 01406 01407 IF (ll_notfound) WRITE(nulprt,*)'psmile_def_metdata: No corner longitudes found' & 01408 ,'for grid', cga_clim_locator(il_port) 01409 ! 01410 !rv,sgi DEALLOCATE(lon) 01411 ! 01412 !Latitude corners 01413 ! 01414 clstrg = cga_clim_locator(il_port)//crn_clim_latsuf 01415 #ifdef __VERBOSE 01416 WRITE(nulprt,*) 'clstrg = ', clstrg 01417 #endif 01418 ll_notfound = .TRUE. 01419 01420 DO il_i=1,il_nvar 01421 01422 CALL mpp_get_atts(il_vars(il_i),name=cl_varname,siz=lens_lat) 01423 #ifdef __VERBOSE 01424 WRITE(nulprt,*) 'cl_varname = ', cl_varname 01425 #endif 01426 IF(TRIM(cl_varname).EQ.TRIM(clstrg) ) THEN 01427 ll_notfound = .FALSE. 01428 #ifdef __VERBOSE 01429 write(nulprt,*) & 01430 'psmile_def_metada, lat, clstrg: ',clstrg,lens_lat 01431 #endif 01432 nbr_corners(il_port)=lens_lat(3) 01433 lla1_tmp = lens_lat(1) 01434 lla2_tmp = lens_lat(2) 01435 ALLOCATE(rl_latscrn(lens_lat(1),lens_lat(2) & 01436 ,nbr_corners(il_port))) 01437 ALLOCATE(latcrn(nbr_corners(il_port) & 01438 ,1:il_xend-il_xbegin+1 & 01439 ,1:il_yend-il_ybegin+1)) 01440 lens_lat(1)=il_xend-il_xbegin+1 01441 lens_lat(2)=il_yend-il_ybegin+1 01442 lens_lat(3)=nbr_corners(il_port) 01443 01444 CALL mpp_read(nc_grdid,il_vars(il_i),rl_latscrn,tindex=1) 01445 !rv,sgi< 01446 if((ig_def_invert(il_port).gt.1).or. & 01447 (ig_def_reverse(il_port).gt.1)) then 01448 llxlon=ig_def_invert(il_port).eq.3 & 01449 .or.ig_def_invert(il_port).eq.4 & 01450 .or.ig_def_reverse(il_port).eq.3 & 01451 .or.ig_def_reverse(il_port).eq.4 01452 llylat=ig_def_invert(il_port).eq.2 & 01453 .or.ig_def_invert(il_port).eq.4 & 01454 .or.ig_def_reverse(il_port).eq.2 & 01455 .or.ig_def_reverse(il_port).eq.4 01456 call inv2dcrn (rl_latscrn,lla1_tmp,lla2_tmp & 01457 ,nbr_corners(il_port),llxlon,llylat) 01458 endif 01459 !rv,sgi> 01460 01461 call transp_crn2cf(rl_latscrn,latcrn,lens_lat) 01462 call mpp_write(ig_units(il_port),crnlatij(il_port),latcrn) 01463 DEALLOCATE(rl_latscrn) 01464 DEALLOCATE(latcrn) 01465 ENDIF 01466 ENDDO 01467 01468 IF (ll_notfound) WRITE(nulprt,*)'psmile_def_metdata: No corner latitudes found' & 01469 ,'for grid', cga_clim_locator(il_port) 01470 01471 DEALLOCATE(lat) 01472 !rv,sgi< 01473 if(ALLOCATED(lon)) DEALLOCATE(lon) 01474 !rv,sgi> 01475 DEALLOCATE(il_vars) 01476 ENDIF 01477 01478 CALL mpp_flush(ig_units(il_port)) 01479 ENDIF 01480 ENDIF 01481 ENDDO 01482 deallocate(c_axis) 01483 deallocate(nbr_corners) 01484 deallocate(crnlonij) 01485 deallocate(crnlatij) 01486 #ifdef __VERBOSE 01487 write(nulprt,*)'| | | Leaving -- psmile_def_metada' 01488 write(nulprt,*)' ' 01489 call flush(nulprt) 01490 #endif 01491 end subroutine psmile_def_metadata 01492 01493 subroutine psmile_close_files(id_error) 01494 !------------------------------------------------------------------------------ 01495 ! This subroutine closes all opened psmile_io files 01496 ! Supposed to be called in prism_terminate before psmile_io_cleanup 01497 ! and before shutting down MPI! 01498 !----------------------------------------------------------------------- 01499 ! 01500 use mod_prism_proto 01501 use mod_comprism_proto 01502 use mod_psmile_io 01503 implicit none 01504 01505 #include <mpif.h> 01506 integer(kind=ip_intwp_p),intent(out)::id_error 01507 !Local variables 01508 integer(kind=ip_intwp_p)::il_port 01509 01510 if( ig_color.ne.1000*ig_mynummod) return 01511 #ifdef __VERBOSE 01512 write(nulprt,*)' ' 01513 write(nulprt,*)'| | | Entering: psmile_close_files' 01514 call flush(nulprt) 01515 #endif 01516 id_error=CLIM_Ok 01517 01518 do il_port=1,nports 01519 ! IF (ig_def_state(il_port) .eq. ip_ignout .or. & 01520 ! ig_def_state(il_port) .eq. ip_expout .or. & 01521 ! ig_def_state(il_port) .eq. ip_output .or. & 01522 ! ig_def_state(il_port) .eq. ip_input ) THEN 01523 !RV:22.01.2003 01524 !Just in case that ig_def_state has been deallocated during termination. 01525 !This way has less interference with rest of prism_terminate. 01526 IF (ig_units(il_port).gt.0) THEN 01527 IF (allocated(dg_times).and. & 01528 associated(dg_times(il_port)%rl_times)) THEN 01529 01530 deallocate(dg_times(il_port)%rl_times) 01531 01532 ENDIF 01533 01534 call mpp_close(ig_units(il_port)) 01535 01536 ENDIF 01537 enddo 01538 01539 !RV: 18.12.2002: Moved these lines into psmile_io_cleanup 01540 ! call mpp_io_exit() 01541 ! call mpp_domains_exit() 01542 ! call mpp_exit() 01543 01544 #ifdef __VERBOSE 01545 WRITE(nulprt,*)'| | | Leaving -- : psmile_close_files' 01546 write(nulprt,*)' ' 01547 call flush(nulprt) 01548 #endif 01549 01550 end subroutine psmile_close_files 01551 01552 subroutine psmile_io_cleanup(id_error) 01553 !----------------------------------------------------------------------- 01554 ! Supposed to be called in prism_terminate after psmile_close_files 01555 !----------------------------------------------------------------------- 01556 use mod_prism_proto 01557 use mod_comprism_proto 01558 use mod_psmile_io 01559 01560 implicit none 01561 01562 integer(kind=ip_intwp_p),intent(out)::id_error 01563 #ifdef __VERBOSE 01564 write(nulprt,*)' ' 01565 write(nulprt,*)'| | | Entering: psmile_io_cleanup' 01566 call flush(nulprt) 01567 #endif 01568 if( ig_color.ne.1000*ig_mynummod) return 01569 !Global arrays to specify IO domains 01570 if(allocated(domain))deallocate(domain) 01571 if(allocated(xdom))deallocate(xdom) 01572 if(allocated(ydom))deallocate(ydom) 01573 if(allocated(lg_is_odd_distribution))deallocate(lg_is_odd_distribution) 01574 01575 !Global arrarys for opening files) 01576 if(allocated(ig_units))deallocate(ig_units) 01577 if(allocated(ig_action))deallocate(ig_action) 01578 if(allocated(ig_form))deallocate(ig_form) 01579 if(allocated(ig_threading))deallocate(ig_threading) 01580 if(allocated(ig_fileset))deallocate(ig_fileset) 01581 01582 !Global arrays for netcdf metadata informations) 01583 if(allocated(ig_var_index))deallocate(ig_var_index) 01584 if(allocated(ig_ntimes))deallocate(ig_ntimes) 01585 if(allocated(dg_times))deallocate( dg_times) 01586 if(allocated(ig_vars))deallocate(ig_vars) 01587 if(allocated(field))deallocate(field) 01588 if(allocated(latij))deallocate(latij) 01589 if(allocated(lonij))deallocate(lonij) 01590 if(allocated(x_axis))deallocate(x_axis) 01591 if(allocated(y_axis))deallocate(y_axis) 01592 if(allocated(t_axis))deallocate(t_axis) 01593 01594 !Was not allocated. Oasis sends 2D arrays. In case of bundles or 3d arrays 01595 !we have to declare z_axis + plus one for bundles. 01596 if(allocated(z_axis))deallocate(z_axis) 01597 01598 ! These arrays are supposed to be part of a different module!) 01599 ! However, no information is available in Oasis 3.0. So, let them live here) 01600 ! for a while.) 01601 01602 if(allocated(cports_lgname)) deallocate(cports_lgname) 01603 if(allocated(cports_units)) deallocate(cports_units) 01604 01605 call mpp_io_exit() 01606 call mpp_domains_exit() 01607 call mpp_exit() 01608 01609 #ifdef __VERBOSE 01610 write(nulprt,*)'| | | Leaving : psmile_io_cleanup' 01611 write(nulprt,*)' ' 01612 call flush(nulprt) 01613 #endif 01614 end subroutine psmile_io_cleanup 01615 01616 subroutine psmile_read_8(id_port_id,rd_field,id_newtime) 01617 !-------------------------------------------------------------------------- 01618 ! Called in rism_get 01619 use mod_kinds_model 01620 use mod_prism_proto 01621 use mod_comprism_proto 01622 use mod_psmile_io 01623 01624 implicit none 01625 01626 integer(kind=ip_intwp_p),intent(in)::id_newtime,id_port_id 01627 REAL(kind=ip_double_p),DIMENSION(myport(4,id_port_id)), 01628 intent(out) :: rd_field 01629 !local 01630 real(kind=ip_double_p)::rl_time 01631 integer(kind=ip_intwp_p)::il_record_no 01632 real(kind=ip_double_p)::diff_first,diff_i 01633 integer(kind=ip_intwp_p)::il_i,i 01634 01635 #ifdef __VERBOSE 01636 write(nulprt,*)' ' 01637 WRITE(nulprt,*) '| | | Entering --- psmile_read_8' 01638 call flush(nulprt) 01639 #endif 01640 01641 !RV:22.01.2003 01642 01643 #ifdef __READ_BY_LOOKUP 01644 !I try to find the record number of timestamp id_newtime by scrolling through 01645 !the time stamps of the file. 01646 01647 rl_time=id_newtime 01648 if(.not.associated(dg_times(id_port_id)%rl_times))then 01649 01650 allocate(dg_times(id_port_id)%rl_times(1:ig_ntimes(id_port_id))) 01651 call mpp_get_times( ig_units(id_port_id) & 01652 ,dg_times(id_port_id)%rl_times(:)) 01653 endif 01654 ! il_record_no=minloc(abs(dg_times(id_port_id)%rl_times-rl_time),dim=1) 01655 il_record_no=1 01656 diff_first=abs(dg_times(id_port_id)%rl_times(1)-rl_time) 01657 do i=1,ig_ntimes(id_port_id) 01658 diff_i=abs(dg_times(id_port_id)%rl_times(i)-rl_time) 01659 if(abs(diff_i).lt.diff_first) then 01660 diff_first=diff_i 01661 il_record_no=i 01662 endif 01663 enddo 01664 if(abs(dg_times(id_port_id)%rl_times(il_record_no)-rl_time).gt.1.d-8) & 01665 then 01666 write(nulprt,*)'psmile_read_8: The time stamp ',id_newtime, & 01667 ' is not available in file ',cg_def_inpfile(id_port_id) 01668 call flush(nulprt) 01669 call MPI_ABORT(mpi_comm,0,mpi_err) 01670 endif 01671 01672 #ifdef __VERBOSE 01673 write(nulprt,*)'psmile_read_8: timestamp, record no.',id_newtime, & 01674 il_record_no 01675 call flush(nulprt) 01676 #endif 01677 #else 01678 !Alternative: 01679 !Assuming that time stamps are starting from zero 01680 il_record_no=id_newtime/ig_def_freq(id_port_id) +1 01681 #endif 01682 01683 call mpp_read(ig_units(id_port_id),ig_vars(id_port_id) & 01684 ,domain(id_port_id),rd_field,il_record_no) 01685 01686 #ifdef __VERBOSE 01687 write(nulprt,*) '| | | Leaving psmile_read_8: End' 01688 write(nulprt,*)' ' 01689 call flush(nulprt) 01690 #endif 01691 end subroutine psmile_read_8 01692 01693 subroutine psmile_read_4(id_port_id,rd_field,id_newtime) 01694 !-------------------------------------------------------------------------- 01695 ! Called in prism_get 01696 use mod_kinds_model 01697 use mod_prism_proto 01698 use mod_comprism_proto 01699 use mod_psmile_io 01700 01701 implicit none 01702 01703 integer(kind=ip_intwp_p),intent(in)::id_newtime,id_port_id 01704 REAL(kind=ip_single_p), DIMENSION(myport(4,id_port_id)) :: rd_field 01705 !local 01706 REAL(kind=ip_double_p), DIMENSION(myport(4,id_port_id)) :: rl_field 01707 real(kind=ip_double_p)::rl_time 01708 integer(kind=ip_intwp_p)::il_record_no,il_i 01709 01710 #ifdef __VERBOSE 01711 write(nulprt,*)' ' 01712 write(nulprt,*) '| | | Entering psmile_read_4 ' 01713 call flush(nulprt) 01714 #endif 01715 01716 !RV:22.01.2003 01717 #ifdef __READ_BY_LOOKUP 01718 !I try to find the record number of timestamp id_newtime by scrolling through 01719 !the time stamps of the file. 01720 rl_time=id_newtime 01721 01722 if(.not.associated(dg_times(id_port_id)%rl_times))then 01723 01724 allocate(dg_times(id_port_id)%rl_times(1:ig_ntimes(id_port_id))) 01725 call mpp_get_times( ig_units(id_port_id), & 01726 dg_times(id_port_id)%rl_times(:)) 01727 01728 endif 01729 01730 il_record_no=minloc(abs(dg_times(id_port_id)%rl_times-rl_time),dim=1) 01731 01732 if(abs(dg_times(id_port_id)%rl_times(il_record_no)-rl_time).gt.1.d-8) & 01733 then 01734 write(nulprt,*)'psmile_read_4: The time stamp ',id_newtime, & 01735 ' is not available in file ',cg_def_inpfile(id_port_id) 01736 call flush(nulprt) 01737 call MPI_ABORT(mpi_comm,0,mpi_err) 01738 endif 01739 #ifdef __VERBOSE 01740 write(nulprt,*)'psmile_read_4: timestamp, record no.',id_newtime, & 01741 il_record_no 01742 call flush(nulprt) 01743 #endif 01744 #else 01745 !Alternative: 01746 !Assuming that time stamps are starting from zero 01747 il_record_no=id_newtime/ig_def_freq(id_port_id) +1 01748 #endif 01749 01750 call mpp_read(ig_units(id_port_id),ig_vars(id_port_id) & 01751 ,domain(id_port_id),rl_field,il_record_no) 01752 rd_field=rl_field 01753 01754 #ifdef __VERBOSE 01755 write(nulprt,*) '| | | Leaving psmile_read_4: End' 01756 write(nulprt,*)' ' 01757 call flush(nulprt) 01758 #endif 01759 end subroutine psmile_read_4 01760 01761 subroutine psmile_write_4(id_port_id,rd_field,id_newtime) 01762 !-------------------------------------------------------------------------- 01763 ! Called in prism_get or prism_put 01764 use mod_kinds_model 01765 use mod_prism_proto 01766 use mod_comprism_proto 01767 use mod_psmile_io 01768 01769 implicit none 01770 integer(kind=ip_intwp_p),intent(in)::id_newtime,id_port_id 01771 REAL(kind=ip_single_p), DIMENSION(myport(4,id_port_id)) :: rd_field 01772 !Local 01773 real(kind=ip_double_p)::rl_time 01774 real(ip_double_p),DIMENSION(myport(4,id_port_id)) ::rl_field 01775 01776 #ifdef __VERBOSE 01777 write(nulprt,*)' ' 01778 write(nulprt,*)'| | | Entering psmile_write_4 ' 01779 call flush(nulprt) 01780 #endif 01781 rl_time=id_newtime 01782 rl_field=rd_field 01783 call mpp_write(ig_units(id_port_id),field(id_port_id) & 01784 ,domain(id_port_id),rl_field,rl_time) 01785 #ifdef __VERBOSE 01786 write(nulprt,*)'| | | Leaving psmile_write_4: End' 01787 write(nulprt,*)' ' 01788 call flush(nulprt) 01789 #endif 01790 01791 end subroutine psmile_write_4 01792 01793 subroutine psmile_write_8(id_port_id,rd_field,id_newtime) 01794 !-------------------------------------------------------------------------- 01795 ! Called in prism_get or prism_put 01796 use mod_kinds_model 01797 use mod_prism_proto 01798 use mod_comprism_proto 01799 use mod_psmile_io 01800 01801 implicit none 01802 integer(kind=ip_intwp_p),intent(in)::id_newtime,id_port_id 01803 integer(kind=ip_intwp_p)::il_stat 01804 REAL(kind=ip_double_p), 01805 DIMENSION(1:myport(4,id_port_id)),intent(inout) :: rd_field 01806 !Local 01807 real(kind=ip_double_p)::rl_time 01808 integer(kind=ip_intwp_p) :: ij 01809 01810 ! write(70+mpi_rank+(ig_color/1000-1)*4,'(4e12.5)')rd_field 01811 #ifdef __VERBOSE 01812 write(nulprt,*)' ' 01813 write(nulprt,*)' | | | Entering psmile_write_8' 01814 call flush(nulprt) 01815 #endif 01816 rl_time=id_newtime 01817 #ifdef __VERBOSE 01818 write(nulprt,*)'rl_time ',rl_time 01819 #endif 01820 ! do ij =1,4608,4608/10 01821 ! write(nulprt,*)'rd_field bef mpp_write ',rd_field(ij) 01822 ! enddo 01823 ! call flush(nulprt) 01824 call mpp_write(ig_units(id_port_id),field(id_port_id) & 01825 ,domain(id_port_id),rd_field,rl_time) 01826 01827 #ifdef __VERBOSE 01828 write(nulprt,*)'| | | Leaving psmile_write_8: End' 01829 write(nulprt,*)' ' 01830 call flush(nulprt) 01831 #endif 01832 end subroutine psmile_write_8 01833 01834 subroutine indexi(n,arr,indx) 01835 !--------------------------------------------------------------------- 01836 ! Generates a list indx which sorts arr in ascending order. 01837 ! This is called sorting by indexing. 01838 ! Code taken from 'Numerical Recipes' 01839 !--------------------------------------------------------------------- 01840 use mod_kinds_model 01841 INTEGER(KIND=IP_INTWP_P),intent(in):: n 01842 INTEGER(KIND=IP_INTWP_P),intent(out)::indx(n) 01843 integer(kind=ip_intwp_p),intent(in):: arr(n) 01844 !Local 01845 integer(kind=ip_intwp_p),PARAMETER:: M=7,NSTACK=128 01846 INTEGER(KIND=IP_INTWP_P):: i,indxt,ir,itemp,j,jstack,k,l,istack(NSTACK) 01847 integer(kind=ip_intwp_p):: a 01848 01849 do 11 j=1,n 01850 indx(j)=j 01851 11 continue 01852 jstack=0 01853 l=1 01854 ir=n 01855 1 if(ir-l.lt.M)then 01856 do 13 j=l+1,ir 01857 indxt=indx(j) 01858 a=arr(indxt) 01859 do 12 i=j-1,1,-1 01860 if(arr(indx(i)).le.a)goto 2 01861 indx(i+1)=indx(i) 01862 12 continue 01863 i=0 01864 2 indx(i+1)=indxt 01865 13 continue 01866 if(jstack.eq.0)return 01867 ir=istack(jstack) 01868 l=istack(jstack-1) 01869 jstack=jstack-2 01870 else 01871 k=(l+ir)/2 01872 itemp=indx(k) 01873 indx(k)=indx(l+1) 01874 indx(l+1)=itemp 01875 if(arr(indx(l+1)).gt.arr(indx(ir)))then 01876 itemp=indx(l+1) 01877 indx(l+1)=indx(ir) 01878 indx(ir)=itemp 01879 endif 01880 if(arr(indx(l)).gt.arr(indx(ir)))then 01881 itemp=indx(l) 01882 indx(l)=indx(ir) 01883 indx(ir)=itemp 01884 endif 01885 if(arr(indx(l+1)).gt.arr(indx(l)))then 01886 itemp=indx(l+1) 01887 indx(l+1)=indx(l) 01888 indx(l)=itemp 01889 endif 01890 i=l+1 01891 j=ir 01892 indxt=indx(l) 01893 a=arr(indxt) 01894 3 continue 01895 i=i+1 01896 if(arr(indx(i)).lt.a)goto 3 01897 4 continue 01898 j=j-1 01899 if(arr(indx(j)).gt.a)goto 4 01900 if(j.lt.i)goto 5 01901 itemp=indx(i) 01902 indx(i)=indx(j) 01903 indx(j)=itemp 01904 goto 3 01905 5 indx(l)=indx(j) 01906 indx(j)=indxt 01907 jstack=jstack+2 01908 if(jstack.gt.NSTACK)pause 'NSTACK too small in indexx' 01909 if(ir-i+1.ge.j-l)then 01910 istack(jstack)=ir 01911 istack(jstack-1)=i 01912 ir=j-1 01913 else 01914 istack(jstack)=j-1 01915 istack(jstack-1)=l 01916 l=i 01917 endif 01918 endif 01919 goto 1 01920 END subroutine indexi 01921 01922 subroutine combine_with_date(cd_in,cd_mode,id_initial_date,cd_on) 01923 !------------------------------------------------------------- 01924 use mod_kinds_model 01925 character(len=*),intent(in)::cd_in 01926 character(len=*),intent(in)::cd_mode 01927 character(len=*),intent(out)::cd_on 01928 integer(kind=ip_intwp_p),intent(in)::id_initial_date(:) 01929 !Local 01930 integer(kind=ip_intwp_p),parameter::il_nsuffixes=2 01931 integer(kind=ip_intwp_p)::leni,leno,lens,lenb,il_i 01932 character(len=8)::suffixes(il_nsuffixes),suffix 01933 character(len=len(cd_in))::cl_basename 01934 character(len=4)::cl_year 01935 character(len=2)::cl_date(1:size(id_initial_date)-1) 01936 01937 suffixes(1)='.nc' 01938 suffixes(2)='.grib' 01939 01940 leni=len(trim(cd_in)) 01941 leno=len(cd_on) 01942 01943 do il_i=1,il_nsuffixes 01944 01945 suffix=suffixes(il_i) 01946 lens=len(trim(suffix)) 01947 01948 if(cd_in(leni-lens+1:leni).eq.suffix(1:lens)) then 01949 cl_basename(1:leni-lens)=cd_in(1:leni-lens) 01950 lenb=leni-lens 01951 EXIT 01952 else 01953 cl_basename(1:leni)=cd_in(1:leni) 01954 lenb=leni 01955 suffix='.nc' 01956 endif 01957 01958 enddo 01959 01960 01961 write(cl_year,'(i4.4)')id_initial_date(1) 01962 01963 do il_i=2,size(id_initial_date) 01964 01965 write(cl_date(il_i-1),'(i2.2)')id_initial_date(il_i) 01966 01967 enddo 01968 01969 if(trim(cd_mode) .eq. 'since' ) then 01970 cd_on=trim(cl_basename(1:lenb))//' '//trim(cd_mode)//' ' & 01971 //cl_year//'-'//cl_date(1)//'-' & 01972 //cl_date(2)//' ' & 01973 //cl_date(3)//':' & 01974 //cl_date(4)//':' & 01975 //cl_date(5) 01976 else 01977 cd_on=trim(cl_basename(1:lenb))//'_'//trim(cd_mode)//'.' & 01978 //cl_year//'-'//cl_date(1)//'-' & 01979 //cl_date(2)//'T' & 01980 //cl_date(3)//':' & 01981 //cl_date(4)//':' & 01982 //cl_date(5)//suffix 01983 endif 01984 01985 return 01986 end subroutine combine_with_date 01987 01988 subroutine transp_crn2cf(x,y,vshape) 01989 use mod_kinds_model 01990 implicit none 01991 ! Input parameters 01992 Integer, Intent(in):: vshape(3) 01993 REAL(kind=ip_double_p),Intent(in)::x(vshape(1),vshape(2),vshape(3)) 01994 ! Ouput parameters 01995 REAL(kind=ip_double_p),Intent(out)::y(vshape(3),vshape(1),vshape(2)) 01996 !Local 01997 integer::il_i,il_j,il_k 01998 01999 do il_k=1,vshape(3) 02000 do il_j=1,vshape(2) 02001 do il_i=1,vshape(1) 02002 y(il_k, il_i, il_j)=x(il_i ,il_j, il_k) 02003 enddo 02004 enddo 02005 enddo 02006 02007 end subroutine transp_crn2cf 02008 02009 !rv,sgi< 02010 subroutine inv2d(x,ilon,ilat,lxlon,lylat) 02011 use mod_kinds_model 02012 use mod_comprism_proto 02013 implicit none 02014 ! 02015 ! Input Parameters 02016 ! 02017 Integer,Intent(In)::ilon,ilat 02018 REAL(kind=ip_double_p)::x(ilon,ilat) 02019 Logical,Intent(In)::lxlon,lylat 02020 ! 02021 ! Local Variables 02022 ! 02023 REAL(kind=ip_double_p)::xtmp 02024 Integer(kind=ip_double_p):: ijmed,iimed,ji,jj 02025 #ifdef __VERBOSE 02026 write(nulprt,*)' ' 02027 write(nulprt,*)'| | | Entering inv2d Start:',ilon,ilat,lxlon,lylat 02028 call flush(nulprt) 02029 #endif 02030 ! 02031 ! Reorder lattitude related values stored as the y-direction. 02032 ! 02033 if(lylat) then 02034 ijmed = ilat/2 02035 DO jj = 1, ijmed 02036 DO ji = 1, ilon 02037 xtmp = x(ji,ilat + 1 - jj) 02038 x(ji,ilat + 1 - jj) = x(ji,jj) 02039 x(ji,jj) = xtmp 02040 ENDDO 02041 ENDDO 02042 endif 02043 02044 ! 02045 ! Reorder longitude related values stored as the x-direction. 02046 ! 02047 if(lxlon) then 02048 iimed = ilon/2 02049 DO jj = 1, ilat 02050 DO ji = 1, iimed 02051 xtmp = x(ilon + 1 - ji,jj) 02052 x(ilon + 1 - ji,jj) = x(ji,jj) 02053 x(ji,jj) = xtmp 02054 ENDDO 02055 ENDDO 02056 endif 02057 02058 #ifdef __VERBOSE 02059 write(nulprt,*)'| | | Leaving inv2d end' 02060 write(nulprt,*)' ' 02061 call flush(nulprt) 02062 #endif 02063 02064 end subroutine inv2d 02065 02066 subroutine inv2dcrn(x,ilon,ilat,icrn,lxlon,lylat) 02067 use mod_kinds_model 02068 use mod_comprism_proto 02069 implicit none 02070 ! 02071 ! Input Parameters 02072 ! 02073 Integer,Intent(In)::ilon,ilat,icrn 02074 REAL(kind=ip_double_p)::x(ilon,ilat,icrn) 02075 Logical,Intent(In)::lxlon,lylat 02076 ! 02077 ! Local Variables 02078 ! 02079 REAL(kind=ip_double_p)::xtmp 02080 Integer(kind=ip_double_p):: ijmed,iimed,ji,jj,jc 02081 #ifdef __VERBOSE 02082 write(nulprt,*)' ' 02083 write(nulprt,*)'| | | Entering inv2dcrn :',ilon,ilat,icrn,lxlon,lylat 02084 call flush(nulprt) 02085 #endif 02086 ! 02087 ! Reorder corners of lattitude 02088 ! 02089 if(lylat) then 02090 ijmed = ilat/2 02091 DO jc=1,icrn 02092 DO jj = 1, ijmed 02093 DO ji = 1, ilon 02094 xtmp = x(ji,ilat + 1 - jj,jc) 02095 x(ji,ilat + 1 - jj,jc) = x(ji,jj,jc) 02096 x(ji,jj,jc) = xtmp 02097 ENDDO 02098 ENDDO 02099 ENDDO 02100 endif 02101 02102 ! 02103 ! Reorder corners of longitude 02104 ! 02105 if(lxlon) then 02106 iimed = ilon/2 02107 DO jc=1,icrn 02108 DO jj = 1, ilat 02109 DO ji = 1, iimed 02110 xtmp = x(ilon + 1 - ji,jj,jc) 02111 x(ilon + 1 - ji,jj,jc) = x(ji,jj,jc) 02112 x(ji,jj,jc) = xtmp 02113 ENDDO 02114 ENDDO 02115 ENDDO 02116 endif 02117 02118 02119 #ifdef __VERBOSE 02120 write(nulprt,*)'| | | Leaving inv2dcrn' 02121 call flush(nulprt) 02122 #endif 02123 end subroutine inv2dcrn 02124 !rv,sgi> 02125 #endif