Back to home page

MITgcm

 
 

    


File indexing completed on 2023-07-14 05:10:44 UTC

view on githubraw file Latest commit de57a2ec on 2023-07-13 16:55:13 UTC
24462d2fa8 Patr*0001 #include "PROFILES_OPTIONS.h"
57c22ecc45 Jean*0002 #include "AD_CONFIG.h"
6a770e0a24 Patr*0003 
d5aa75d39a Jean*0004 C     *==========================================================*
                0005 C     | subroutine profiles_init_fixed
                0006 C     | o initialization for netcdf profiles data
                0007 C     | started: Gael Forget 15-March-2006
                0008 C     | extended: Gael Forget 14-June-2007
                0009 C     *==========================================================*
6a770e0a24 Patr*0010 
                0011       SUBROUTINE profiles_init_fixed( myThid )
                0012 
                0013       implicit none
                0014 
                0015 C ==================== Global Variables ===========================
                0016 #include "SIZE.h"
                0017 #include "EEPARAMS.h"
                0018 #include "PARAMS.h"
                0019 #include "GRID.h"
                0020 #include "DYNVARS.h"
d28c90138c Patr*0021 #ifdef ALLOW_CAL
6a770e0a24 Patr*0022 #include "cal.h"
d28c90138c Patr*0023 #endif
24462d2fa8 Patr*0024 #ifdef ALLOW_PROFILES
6328b73337 Gael*0025 # include "PROFILES_SIZE.h"
6e4c90fea3 Patr*0026 # include "profiles.h"
                0027 # include "netcdf.inc"
                0028 #endif
39ce977435 Gael*0029 
                0030 C ==================== Routine Arguments ==========================
                0031 
                0032       integer myThid
                0033 
6a770e0a24 Patr*0034 C ==================== Routine Variables ==========================
                0035 
39ce977435 Gael*0036 #ifdef ALLOW_PROFILES
                0037 
a996b8bdc0 Gael*0038       integer i,j,k,l,m,bi,bj,iG,jG,num_file,ProfNo_tile
c9c84c2afb Gael*0039       integer stopProfiles
6a770e0a24 Patr*0040       integer fid, dimid, varid1, varid1a, varid1b
                0041       integer varid2,varid3
                0042       _RL tmpyymmdd(1000),tmphhmmss(1000),diffsecs
5042c05de8 Gael*0043       _RL yymmddMin,yymmddMax
                0044       _RL hhmmssMin,hhmmssMax
                0045 
b2a948f981 Gael*0046       integer tmpdate(4),tmpdiff(4),profIsInRunTime
39ce977435 Gael*0047       _RL  tmp_lon, tmp_lon2(1000), tmp_lat2(1000), lon_cur, lat_cur
                0048       _RL lon_1, lon_2, lat_1, lat_2
2767dff983 Jean*0049       _RL lon_tmp1, lon_tmp2
39ce977435 Gael*0050       _RL lat_fac, lon_fac
                0051       integer prof_i, prof_j
6a770e0a24 Patr*0052       integer vec_start(2), vec_count(2), profno_div1000, kk
de57a2ec4b Mart*0053       character*(MAX_LEN_FNAM) profilesfile, fnamedatanc
                0054       character*(MAX_LEN_FNAM) fnameequinc, adfnameequinc
38287224dd Gael*0055       integer IL, JL, KL, err
6a770e0a24 Patr*0056       logical  exst
                0057 
ba63501b4c Gael*0058       integer varid_intp1, varid_intp2, varid_intp11 , varid_intp22
38287224dd Gael*0059       integer varid_intp3, varid_intp4, varid_intp5, q, iINTERP
ba63501b4c Gael*0060       _RL tmp_i(1000,NUM_INTERP_POINTS)
                0061       _RL tmp_j(1000,NUM_INTERP_POINTS)
                0062       _RL tmp_weights(1000,NUM_INTERP_POINTS),tmp_sum_weights
                0063       _RL tmp_xC11(1000),tmp_yC11(1000)
                0064       _RL tmp_xCNINJ(1000),tmp_yCNINJ(1000)
c9c84c2afb Gael*0065       integer stopGenericGrid
ba63501b4c Gael*0066       Real*8 xy_buffer_r8(0:sNx+1,0:sNy+1)
                0067       integer vec_start2(2), vec_count2(2)
a996b8bdc0 Gael*0068       integer hh, ProfNo_hh
6b2230d510 Ou W*0069 #ifdef ALLOW_PROFILES_SAMPLESPLIT_COST
                0070       integer varid4
                0071       _RL tmp_avgbin(1000)
                0072 #endif
ba63501b4c Gael*0073 
6a770e0a24 Patr*0074 c     == external functions ==
                0075       integer ILNBLNK
b2a948f981 Gael*0076       EXTERNAL ILNBLNK
ba63501b4c Gael*0077       integer MDS_RECLEN
b2a948f981 Gael*0078       EXTERNAL MDS_RECLEN
ce2e1d3cd5 Patr*0079       character*(max_len_mbuf) msgbuf
6a770e0a24 Patr*0080 
                0081 c--   == end of interface ==
                0082 
f0e4bffe35 Gael*0083       write(msgbuf,'(a)') ' '
                0084       call print_message( msgbuf, standardmessageunit,
                0085      &                    SQUEEZE_RIGHT , mythid)
                0086       write(msgbuf,'(a)')
                0087      &'// ======================================================='
                0088       call print_message( msgbuf, standardmessageunit,
                0089      &                    SQUEEZE_RIGHT , mythid)
                0090       write(msgbuf,'(a)')
                0091      &'// insitu profiles model sampling >>> START <<<'
                0092       call print_message( msgbuf, standardmessageunit,
                0093      &                    SQUEEZE_RIGHT , mythid)
                0094       write(msgbuf,'(a)')
                0095      &'// ======================================================='
                0096       call print_message( msgbuf, standardmessageunit,
                0097      &                    SQUEEZE_RIGHT , mythid)
                0098       write(msgbuf,'(a)') ' '
                0099       call print_message( msgbuf, standardmessageunit,
                0100      &                    SQUEEZE_RIGHT , mythid)
                0101 
c9c84c2afb Gael*0102       stopProfiles=0
                0103       stopGenericGrid=0
f0e4bffe35 Gael*0104 
                0105       IF ( (.NOT.profilesDoGenGrid).AND.
                0106      &     (.NOT.usingSphericalPolarGrid .OR. rotateGrid) ) THEN
                0107         WRITE(msgBuf,'(2A)') 'PROFILES_INIT_FIXED: ',
                0108      &  'profilesDoGenGrid=.true. is required'
                0109         CALL PRINT_ERROR( msgBuf , myThid )
                0110         WRITE(msgBuf,'(2A)') 'PROFILES_INIT_FIXED: ',
                0111      &  'unless usingSphericalGrid=.TRUE. and rotateGrid=.FALSE.'
                0112         CALL PRINT_ERROR( msgBuf , myThid )
b00d6c1700 Gael*0113         CALL ALL_PROC_DIE( myThid )
f0e4bffe35 Gael*0114         STOP 'ABNORMAL END: S/R PROFILES_INIT_FIXED'
                0115       ENDIF
                0116 
                0117       write(msgbuf,'(a)') ' '
                0118       call print_message( msgbuf, standardmessageunit,
                0119      &                    SQUEEZE_RIGHT , mythid)
64278b1175 Jean*0120       write(msgbuf,'(a)') 'general packages parameters :'
f0e4bffe35 Gael*0121       JL  = ILNBLNK( profilesDir )
                0122       if (JL.NE.0) then
80468c4098 Matt*0123         write(msgbuf,'(a,a)') '  profilesDir ',profilesDir(1:JL)
f0e4bffe35 Gael*0124       else
64278b1175 Jean*0125         write(msgbuf,'(a,a)') '  profilesDir ','./'
f0e4bffe35 Gael*0126       endif
                0127       call print_message(
                0128      &   msgbuf, standardmessageunit, SQUEEZE_RIGHT , mythid)
64278b1175 Jean*0129       write(msgbuf,'(a,l5)') '  profilesDoGenGrid  ',profilesDoGenGrid
f0e4bffe35 Gael*0130       call print_message(
                0131      &   msgbuf, standardmessageunit, SQUEEZE_RIGHT , mythid)
64278b1175 Jean*0132       write(msgbuf,'(a,l5)') '  profilesDoNcOutput ',profilesDoNcOutput
f0e4bffe35 Gael*0133       call print_message(
                0134      &   msgbuf, standardmessageunit, SQUEEZE_RIGHT , mythid)
                0135       write(msgbuf,'(a)') ' '
                0136       call print_message( msgbuf, standardmessageunit,
                0137      &                    SQUEEZE_RIGHT , mythid)
                0138 
ba63501b4c Gael*0139       _BEGIN_MASTER( mythid )
b00d6c1700 Gael*0140 
ba63501b4c Gael*0141       DO bj=1,nSy
64278b1175 Jean*0142       DO bi=1,nSx
71a5587721 Gael*0143 
a377358620 Patr*0144         profiles_curfile_buff(bi,bj)=0
5042c05de8 Gael*0145         yymmddMin=modelstartdate(1)
                0146         yymmddMax=modelenddate(1)
                0147         hhmmssMin=modelstartdate(2)
                0148         hhmmssMax=modelenddate(2)
6e4c90fea3 Patr*0149 
a377358620 Patr*0150         do m=1,NLEVELMAX
                0151          do l=1,1000
                0152           do k=1,NVARMAX
c9c84c2afb Gael*0153            profiles_data_buff(m,l,k,bi,bj)=0. _d 0
                0154            profiles_weight_buff(m,l,k,bi,bj)=0. _d 0
a377358620 Patr*0155           enddo
                0156          enddo
                0157         enddo
d5aa75d39a Jean*0158 
a377358620 Patr*0159         do num_file=1,NFILESPROFMAX
6a770e0a24 Patr*0160 
a996b8bdc0 Gael*0161       ProfNo_hh=0
                0162 
38287224dd Gael*0163       profilesfile=' '
6a770e0a24 Patr*0164       IL  = ILNBLNK( profilesfiles(num_file) )
                0165       if (IL.NE.0) then
de57a2ec4b Mart*0166         write(profilesfile,'(1a)')
ce2e1d3cd5 Patr*0167      &     profilesfiles(num_file)(1:IL)
38287224dd Gael*0168         write(msgbuf,'(a)') ' '
                0169         call print_message( msgbuf, standardmessageunit,
                0170      &                    SQUEEZE_RIGHT , mythid)
f0e4bffe35 Gael*0171         write(msgbuf,'(a,i3,a,a)')
a996b8bdc0 Gael*0172      &     'profiles file #',num_file,' is ', profilesfile(1:IL)
ce2e1d3cd5 Patr*0173         call print_message(
                0174      &     msgbuf, standardmessageunit, SQUEEZE_RIGHT , mythid)
6a770e0a24 Patr*0175       endif
                0176 
                0177       IL  = ILNBLNK( profilesfile )
                0178       if (IL.NE.0) then
                0179 
                0180 C===========================================================
ba63501b4c Gael*0181 c open data files and read information
6a770e0a24 Patr*0182 C===========================================================
                0183 
de57a2ec4b Mart*0184       write(fnamedatanc,'(2a)') profilesfile(1:IL),'.nc'
71a5587721 Gael*0185       err = NF_OPEN(fnamedatanc, 0, fiddata(num_file,bi,bj))
6a770e0a24 Patr*0186 
                0187 c1)  read the number of profiles :
71a5587721 Gael*0188       fid=fiddata(num_file,bi,bj)
6a770e0a24 Patr*0189       err = NF_INQ_DIMID(fid,'iPROF', dimid )
71a5587721 Gael*0190       err = NF_INQ_DIMLEN(fid, dimid, ProfNo(num_file,bi,bj) )
6a770e0a24 Patr*0191       err = NF_INQ_DIMID(fid,'iDEPTH', dimid )
                0192       if (err.NE.NF_NOERR) then
ce2e1d3cd5 Patr*0193         err = NF_INQ_DIMID(fid,'Z', dimid )
6a770e0a24 Patr*0194       endif
71a5587721 Gael*0195       err = NF_INQ_DIMLEN(fid, dimid, ProfDepthNo(num_file,bi,bj) )
38287224dd Gael*0196       err = NF_INQ_DIMID(fid,'iINTERP', dimid )
                0197       if (err.EQ.NF_NOERR) then
                0198         err = NF_INQ_DIMLEN(fid, dimid, iINTERP )
                0199       else
                0200         iINTERP=NUM_INTERP_POINTS
                0201       endif
f0e4bffe35 Gael*0202 
a996b8bdc0 Gael*0203       write(msgbuf,'(a,i4,a,i4)')
                0204      &   '  current tile is bi,bj                      =',
                0205      &   bi,',',bj
f0e4bffe35 Gael*0206       call print_message(
                0207      &   msgbuf, standardmessageunit, SQUEEZE_RIGHT , mythid)
a996b8bdc0 Gael*0208       write(msgbuf,'(a,i9)')
                0209      &   '  # of depth levels in file                  =',
                0210      &   ProfDepthNo(num_file,bi,bj)
f0e4bffe35 Gael*0211       call print_message(
                0212      &   msgbuf, standardmessageunit, SQUEEZE_RIGHT , mythid)
64278b1175 Jean*0213       write(msgbuf,'(a,i9)')
a996b8bdc0 Gael*0214      &   '  # of profiles in file                      =',
                0215      &   ProfNo(num_file,bi,bj)
ce2e1d3cd5 Patr*0216       call print_message(
                0217      &   msgbuf, standardmessageunit, SQUEEZE_RIGHT , mythid)
6a770e0a24 Patr*0218 
                0219 c2) read the dates and positions :
b632e3ba1b Gael*0220       err = NF_INQ_VARID(fid,'prof_depth', varid1a )
                0221       if (err.NE.NF_NOERR) then
                0222 c       if no prof_depth is found, then try old variable name:
                0223         err = NF_INQ_VARID(fid,'depth', varid1a )
                0224       endif
                0225       if (err.NE.NF_NOERR) then
                0226 c       if neither is found, then stop
f0e4bffe35 Gael*0227         IL  = ILNBLNK( profilesfile )
                0228         WRITE(msgBuf,'(3A)')
                0229      & 'PROFILES_INIT_FIXED: file ', profilesfile(1:IL),
                0230      & '.nc is not in the pkg/profiles format (no prof_depth etc.)'
                0231         CALL PRINT_ERROR( msgBuf, myThid)
c9c84c2afb Gael*0232         stopProfiles=1
b632e3ba1b Gael*0233       endif
                0234 
71a5587721 Gael*0235       do k=1,ProfDepthNo(num_file,bi,bj)
6a770e0a24 Patr*0236       err = NF_GET_VAR1_DOUBLE(fid,varid1a,k,
71a5587721 Gael*0237      & prof_depth(num_file,k,bi,bj))
6a770e0a24 Patr*0238       enddo
                0239 
                0240       err = NF_INQ_VARID(fid,'prof_YYYYMMDD', varid1a )
                0241       err = NF_INQ_VARID(fid,'prof_HHMMSS', varid1b )
                0242       err = NF_INQ_VARID(fid,'prof_lon', varid2 )
                0243       err = NF_INQ_VARID(fid,'prof_lat', varid3 )
6b2230d510 Ou W*0244 #ifdef ALLOW_PROFILES_SAMPLESPLIT_COST
                0245       err = NF_INQ_VARID(fid,'prof_bin_id_a', varid4 )
                0246 #endif
6a770e0a24 Patr*0247 
ba63501b4c Gael*0248       if (err.NE.NF_NOERR) then
f0e4bffe35 Gael*0249         IL  = ILNBLNK( profilesfile )
                0250         WRITE(msgBuf,'(3A)')
                0251      & 'PROFILES_INIT_FIXED: file ', profilesfile(1:IL),
                0252      & '.nc is not in the pkg/profiles format (no prof_YYYYMMDD etc.)'
                0253         CALL PRINT_ERROR( msgBuf, myThid)
c9c84c2afb Gael*0254       stopProfiles=1
ba63501b4c Gael*0255       endif
                0256 
f0e4bffe35 Gael*0257       if (profilesDoGenGrid) then
d5aa75d39a Jean*0258 c3) read interpolattion information (grid points, coeffs, etc.)
ba63501b4c Gael*0259            err = NF_INQ_VARID(fid,'prof_interp_XC11',varid_intp1)
                0260            err = NF_INQ_VARID(fid,'prof_interp_YC11',varid_intp2)
                0261            err = NF_INQ_VARID(fid,'prof_interp_XCNINJ',varid_intp11)
                0262            err = NF_INQ_VARID(fid,'prof_interp_YCNINJ',varid_intp22)
                0263            err = NF_INQ_VARID(fid,'prof_interp_weights',varid_intp3)
                0264            err = NF_INQ_VARID(fid,'prof_interp_i',varid_intp4)
                0265            err = NF_INQ_VARID(fid,'prof_interp_j',varid_intp5)
                0266       if (err.NE.NF_NOERR) then
f0e4bffe35 Gael*0267         IL  = ILNBLNK( profilesfile )
                0268         WRITE(msgBuf,'(3A)')
                0269      & 'PROFILES_INIT_FIXED: file ', profilesfile(1:IL),
                0270      & '.nc is missing interpolation information (profilesDoGenGrid)'
                0271         CALL PRINT_ERROR( msgBuf, myThid)
c9c84c2afb Gael*0272       stopGenericGrid=2
ba63501b4c Gael*0273       endif
f0e4bffe35 Gael*0274       endif
ba63501b4c Gael*0275 
                0276 c4) default values
6a770e0a24 Patr*0277       do k=1,NOBSGLOB
c9c84c2afb Gael*0278       prof_time(num_file,k,bi,bj)=-999. _d 0
                0279       prof_lon(num_file,k,bi,bj)=-999. _d 0
                0280       prof_lat(num_file,k,bi,bj)=-999. _d 0
                0281       prof_ind_glob(num_file,k,bi,bj)=0
6b2230d510 Ou W*0282 #ifdef ALLOW_PROFILES_SAMPLESPLIT_COST
                0283       prof_ind_avgbin(num_file,k,bi,bj)=-999
                0284 #endif
ba63501b4c Gael*0285       do q = 1,NUM_INTERP_POINTS
38287224dd Gael*0286          prof_interp_i(num_file,k,q,bi,bj) = 1
                0287          prof_interp_j(num_file,k,q,bi,bj) = 1
                0288          prof_interp_weights(num_file,k,q,bi,bj) = 0. _d 0
ba63501b4c Gael*0289       enddo
c9c84c2afb Gael*0290       prof_interp_xC11(num_file,k,bi,bj)=-999. _d 0
                0291       prof_interp_yC11(num_file,k,bi,bj)=-999. _d 0
                0292       prof_interp_xCNINJ(num_file,k,bi,bj)=-999. _d 0
                0293       prof_interp_yCNINJ(num_file,k,bi,bj)=-999. _d 0
6a770e0a24 Patr*0294       enddo
                0295 
ba63501b4c Gael*0296 c5) main loop: look for profiles in this tile
a996b8bdc0 Gael*0297       ProfNo_tile=0
71a5587721 Gael*0298       profno_div1000=max(0,int(ProfNo(num_file,bi,bj)/1000))
6a770e0a24 Patr*0299 
f527c11034 Gael*0300       do kk=1,profno_div1000+1
6a770e0a24 Patr*0301 
71a5587721 Gael*0302       if (min(ProfNo(num_file,bi,bj), 1000*kk).GE.
6a770e0a24 Patr*0303      &  1+1000*(kk-1)) then
                0304 
ba63501b4c Gael*0305 c5.1) read a chunk
6a770e0a24 Patr*0306       vec_start(1)=1
                0307       vec_start(2)=1+1000*(kk-1)
                0308       vec_count(1)=1
71a5587721 Gael*0309       vec_count(2)=min(1000,ProfNo(num_file,bi,bj)-1000*(kk-1))
6a770e0a24 Patr*0310 
                0311       if ( (vec_count(2).LE.0).OR.(vec_count(2).GT.1000).OR.
                0312      & (vec_start(2).LE.0).OR.
d5aa75d39a Jean*0313      & (vec_count(2)+vec_start(2)-1.GT.ProfNo(num_file,bi,bj)) )
6a770e0a24 Patr*0314      & then
f0e4bffe35 Gael*0315         IL  = ILNBLNK( profilesfile )
                0316         WRITE(msgBuf,'(3A)')
                0317      & 'PROFILES_INIT_FIXED: file ', profilesfile(1:IL),
                0318      & '.nc was not read properly (case 1).'
                0319         CALL PRINT_ERROR( msgBuf, myThid)
c9c84c2afb Gael*0320       stopProfiles=1
6a770e0a24 Patr*0321       endif
                0322 
                0323       err = NF_GET_VARA_DOUBLE(fid,varid1a,vec_start(2),
                0324      & vec_count(2), tmpyymmdd)
                0325       err = NF_GET_VARA_DOUBLE(fid,varid1b,vec_start(2),
                0326      & vec_count(2), tmphhmmss)
                0327       err = NF_GET_VARA_DOUBLE(fid,varid2,vec_start(2),
                0328      & vec_count(2), tmp_lon2)
                0329       err = NF_GET_VARA_DOUBLE(fid,varid3,vec_start(2),
                0330      & vec_count(2), tmp_lat2)
6b2230d510 Ou W*0331 #ifdef ALLOW_PROFILES_SAMPLESPLIT_COST
                0332       err = NF_GET_VARA_DOUBLE(fid,varid4,vec_start(2),
                0333      & vec_count(2), tmp_avgbin)
                0334 #endif
6a770e0a24 Patr*0335 
                0336       if (err.NE.NF_NOERR) then
f0e4bffe35 Gael*0337         WRITE(msgBuf,'(3A)')
                0338      & 'PROFILES_INIT_FIXED: file ', profilesfile(1:IL),
                0339      & '.nc was not read properly (case 2).'
                0340         CALL PRINT_ERROR( msgBuf, myThid)
c9c84c2afb Gael*0341       stopProfiles=1
d5aa75d39a Jean*0342       endif
6a770e0a24 Patr*0343 
39ce977435 Gael*0344 c if profilesDoGenGrid then also read in the interpolation coeffs and indices
f0e4bffe35 Gael*0345       if (profilesDoGenGrid) then
d5aa75d39a Jean*0346       err = NF_GET_VARA_DOUBLE(fid,varid_intp1,vec_start(2),
ba63501b4c Gael*0347      & vec_count(2), tmp_xC11)
d5aa75d39a Jean*0348       err = NF_GET_VARA_DOUBLE(fid,varid_intp2,vec_start(2),
ba63501b4c Gael*0349      & vec_count(2), tmp_yC11)
                0350       err = NF_GET_VARA_DOUBLE(fid,varid_intp11,vec_start(2),
                0351      & vec_count(2), tmp_xCNINJ)
d5aa75d39a Jean*0352       err = NF_GET_VARA_DOUBLE(fid,varid_intp22,vec_start(2),
ba63501b4c Gael*0353      & vec_count(2), tmp_yCNINJ)
38287224dd Gael*0354       do q=1,iINTERP
ba63501b4c Gael*0355         vec_start2(1)=q
                0356         vec_start2(2)=1+1000*(kk-1)
                0357         vec_count2(1)=1
                0358         vec_count2(2)=min(1000,ProfNo(num_file,bi,bj)-1000*(kk-1))
d5aa75d39a Jean*0359         err = NF_GET_VARA_DOUBLE(fid,varid_intp3,vec_start2,
ba63501b4c Gael*0360      &  vec_count2, tmp_weights(1,q))
d5aa75d39a Jean*0361         err = NF_GET_VARA_DOUBLE(fid,varid_intp4,vec_start2,
ba63501b4c Gael*0362      &  vec_count2, tmp_i(1,q))
d5aa75d39a Jean*0363         err = NF_GET_VARA_DOUBLE(fid,varid_intp5,vec_start2,
ba63501b4c Gael*0364      &  vec_count2, tmp_j(1,q))
                0365       enddo
f0e4bffe35 Gael*0366       endif
ba63501b4c Gael*0367 
                0368 c5.2) loop through this chunk
71a5587721 Gael*0369       do k=1,min(1000,ProfNo(num_file,bi,bj)-1000*(kk-1))
6a770e0a24 Patr*0370 
c9c84c2afb Gael*0371       if ( stopProfiles .EQ. 0) then
ba63501b4c Gael*0372 
b2a948f981 Gael*0373       profIsInRunTime=1
                0374 
a996b8bdc0 Gael*0375       IF (( ( tmpyymmdd(k).GT.yymmddMin ).OR.(( tmpyymmdd(k).EQ.
5042c05de8 Gael*0376      &        yymmddMin ).AND.( tmphhmmss(k).GT.hhmmssMin ))).AND.
                0377      &    ( ( tmpyymmdd(k).LT.yymmddMax ).OR.(( tmpyymmdd(k).EQ.
a996b8bdc0 Gael*0378      &        yymmddMax ).AND.( tmphhmmss(k).LT.hhmmssMax ))) ) THEN
                0379         hh = int(tmphhmmss(k))/10000
                0380         IF ( hh.LT.hoursPerDay ) THEN
                0381           profIsInRunTime=1
                0382           call cal_FullDate( int(tmpyymmdd(k)),int(tmphhmmss(k)),
                0383      &     tmpdate,mythid )
                0384           call cal_TimePassed( modelstartdate,tmpdate,tmpdiff,mythid )
                0385           call cal_ToSeconds (tmpdiff,diffsecs,mythid)
                0386           diffsecs=diffsecs+nIter0*deltaTclock
                0387         ELSE
                0388 c if tmphhmmss is out of range then disregard profile
                0389           profIsInRunTime=0
                0390           diffsecs=-deltaTclock
                0391           ProfNo_hh=ProfNo_hh+1
                0392         ENDIF
                0393       ELSE
b2a948f981 Gael*0394         profIsInRunTime=0
                0395         diffsecs=-deltaTclock
a996b8bdc0 Gael*0396       ENDIF
6a770e0a24 Patr*0397 
39ce977435 Gael*0398 c ==============================================================================
                0399 
                0400 c 5.2a) determine whether profiles is in current tile domain (lat-lon grid case)
                0401        if ((.NOT.profilesDoGenGrid).AND.(profIsInRunTime.EQ.1)) then
                0402 
f527c11034 Gael*0403        if (xC(sNx+1,1,bi,bj).LT.xC(1,1,bi,bj)) then
c9c84c2afb Gael*0404         tmp_lon=xC(sNx+1,1,bi,bj)+360. _d 0
6a770e0a24 Patr*0405        else
f527c11034 Gael*0406         tmp_lon=xC(sNx+1,1,bi,bj)
6a770e0a24 Patr*0407        endif
39ce977435 Gael*0408 
6a770e0a24 Patr*0409        if ((xC(1,1,bi,bj).LE.tmp_lon2(k)).AND.
f527c11034 Gael*0410      & (tmp_lon.GT.tmp_lon2(k)).AND.
                0411      & (yC(1,1,bi,bj).LE.tmp_lat2(k)).AND.
39ce977435 Gael*0412      & (yC(1,sNy+1,bi,bj).GT.tmp_lat2(k)) ) then
                0413          lon_cur=tmp_lon2(k)
                0414          lat_cur=tmp_lat2(k)
                0415        elseif ((xC(sNx+1,1,bi,bj).LT.xC(1,1,bi,bj)).AND.
                0416      &  (xC(1,1,bi,bj).LE.tmp_lon2(k)+360. _d 0).AND.
c9c84c2afb Gael*0417      &  (tmp_lon.GT.tmp_lon2(k)+360. _d 0).AND.
6a770e0a24 Patr*0418      &  (yC(1,1,bi,bj).LE.tmp_lat2(k)).AND.
                0419      &  (yC(1,sNy+1,bi,bj).GT.tmp_lat2(k))
f527c11034 Gael*0420      &  ) then
39ce977435 Gael*0421          lon_cur=tmp_lon2(k)+360. _d 0
                0422          lat_cur=tmp_lat2(k)
f0e4bffe35 Gael*0423        else
39ce977435 Gael*0424          profIsInRunTime=0
                0425        endif
                0426 
                0427 c now determine value of i,j to the south-ouest of data point
                0428        prof_i=-10
                0429        prof_j=-10
                0430        lon_1=-10
                0431        lon_2=-10
                0432        lat_1=-10
                0433        lat_2=-10
                0434 
8da32d4429 Gael*0435        if (profIsInRunTime.EQ.1) then
39ce977435 Gael*0436         DO j=1,sNy+1
                0437          DO i=1,sNx+1
                0438 
                0439 c value of j, south of the data point:
                0440         if ((yC(i,j,bi,bj).LE.lat_cur).AND.
                0441      &      (yC(i,j+1,bi,bj).GT.lat_cur)) then
                0442           prof_j=j
                0443           lat_1=yC(i,j,bi,bj)
                0444           lat_2=yC(i,j+1,bi,bj)
                0445         endif
                0446 
                0447 c value of i, west of the data point:
                0448          if (xC(i+1,j,bi,bj).LT.xC(1,j,bi,bj)) then
                0449            lon_tmp2=xC(i+1,j,bi,bj)+360
                0450          else
                0451            lon_tmp2=xC(i+1,j,bi,bj)
                0452          endif
                0453          if (xC(i,j,bi,bj).LT.xC(1,j,bi,bj)) then
                0454            lon_tmp1=xC(i,j,bi,bj)+360
                0455          else
                0456            lon_tmp1=xC(i,j,bi,bj)
                0457          endif
                0458 
                0459          if ((lon_tmp1.LE.lon_cur).AND.(lon_tmp2.GT.lon_cur)) then
                0460            prof_i=i
                0461            lon_1=lon_tmp1
                0462            lon_2=lon_tmp2
                0463          endif
                0464 
                0465         ENDDO
                0466        ENDDO
8da32d4429 Gael*0467       endif
39ce977435 Gael*0468 
                0469       if ((prof_i.EQ.-10).OR.(prof_j.EQ.-10)) profIsInRunTime=0
                0470 
                0471       if (profIsInRunTime.EQ.1) then
                0472 c if yes then store prof_time and longitude and latitude:
a996b8bdc0 Gael*0473         ProfNo_tile=ProfNo_tile+1
                0474         prof_time(num_file,ProfNo_tile,bi,bj)=diffsecs
                0475         prof_lon(num_file,ProfNo_tile,bi,bj)=lon_cur
                0476         prof_lat(num_file,ProfNo_tile,bi,bj)=lat_cur
                0477         prof_ind_glob(num_file,ProfNo_tile,bi,bj)=k+1000*(kk-1)
6b2230d510 Ou W*0478 #ifdef ALLOW_PROFILES_SAMPLESPLIT_COST
                0479         prof_ind_avgbin(num_file,ProfNo_tile,bi,bj)=tmp_avgbin(k)
                0480 #endif
39ce977435 Gael*0481 c then store interpolation coeffs and indices
                0482         lon_fac=(lon_cur-lon_1)/(lon_2-lon_1)
                0483         lat_fac=(lat_cur-lat_1)/(lat_2-lat_1)
a996b8bdc0 Gael*0484         prof_interp_weights(num_file,ProfNo_tile,1,bi,bj)=
39ce977435 Gael*0485      &     (1-lon_fac)*(1-lat_fac)
a996b8bdc0 Gael*0486         prof_interp_i(num_file,ProfNo_tile,1,bi,bj)=prof_i
                0487         prof_interp_j(num_file,ProfNo_tile,1,bi,bj)=prof_j
                0488         prof_interp_weights(num_file,ProfNo_tile,2,bi,bj)=
39ce977435 Gael*0489      &     lon_fac*(1-lat_fac)
a996b8bdc0 Gael*0490         prof_interp_i(num_file,ProfNo_tile,2,bi,bj)=prof_i+1
                0491         prof_interp_j(num_file,ProfNo_tile,2,bi,bj)=prof_j
                0492         prof_interp_weights(num_file,ProfNo_tile,3,bi,bj)=
39ce977435 Gael*0493      &     (1-lon_fac)*lat_fac
a996b8bdc0 Gael*0494         prof_interp_i(num_file,ProfNo_tile,3,bi,bj)=prof_i
                0495         prof_interp_j(num_file,ProfNo_tile,3,bi,bj)=prof_j+1
                0496         prof_interp_weights(num_file,ProfNo_tile,4,bi,bj)=
39ce977435 Gael*0497      &     lon_fac*lat_fac
a996b8bdc0 Gael*0498         prof_interp_i(num_file,ProfNo_tile,4,bi,bj)=prof_i+1
                0499         prof_interp_j(num_file,ProfNo_tile,4,bi,bj)=prof_j+1
39ce977435 Gael*0500 
                0501       endif
                0502 
                0503 c ==============================================================================
                0504 
                0505 c 5.2a) determine whether profiles is in current tile domain (generic grid case)
                0506 
                0507        elseif (profIsInRunTime.EQ.1) then
                0508 
c9c84c2afb Gael*0509        if (stopGenericGrid.EQ.0) then
ba63501b4c Gael*0510 
c9c84c2afb Gael*0511        if ( ( abs( tmp_xC11(k) - xC(1,1,bi,bj) ).LT.0.0001 _d 0 ) .AND.
                0512      & ( abs( tmp_yC11(k) - yC(1,1,bi,bj) ).LT.0.0001 _d 0) .AND.
                0513      & ( abs( tmp_xCNINJ(k) - xC(sNx,sNy,bi,bj) ).LT.0.0001 _d 0 ) .AND.
                0514      & ( abs( tmp_yCNINJ(k) - yC(sNx,sNy,bi,bj) ).LT.0.0001 _d 0 )
b2a948f981 Gael*0515      & .AND.(profIsInRunTime.EQ.1)) then
ba63501b4c Gael*0516 
39ce977435 Gael*0517 c if yes then store prof_time and interpolation coeffs and indices:
a996b8bdc0 Gael*0518        ProfNo_tile=ProfNo_tile+1
                0519        prof_time(num_file,ProfNo_tile,bi,bj)=diffsecs
6b2230d510 Ou W*0520 #ifdef ALLOW_PROFILES_SAMPLESPLIT_COST
                0521          prof_ind_avgbin(num_file,ProfNo_tile,bi,bj)=tmp_avgbin(k)
                0522 #endif
a996b8bdc0 Gael*0523        prof_interp_xC11(num_file,ProfNo_tile,bi,bj)=tmp_xC11(k)
                0524        prof_interp_yC11(num_file,ProfNo_tile,bi,bj)=tmp_yC11(k)
                0525        prof_interp_xCNINJ(num_file,ProfNo_tile,bi,bj)=tmp_xCNINJ(k)
                0526        prof_interp_yCNINJ(num_file,ProfNo_tile,bi,bj)=tmp_yCNINJ(k)
ba63501b4c Gael*0527        tmp_sum_weights=0. _d 0
38287224dd Gael*0528         do q = 1,iINTERP
a996b8bdc0 Gael*0529              prof_interp_weights(num_file,ProfNo_tile,q,bi,bj)
ba63501b4c Gael*0530      &       =tmp_weights(k,q)
a996b8bdc0 Gael*0531              prof_interp_i(num_file,ProfNo_tile,q,bi,bj)
ba63501b4c Gael*0532      &       =tmp_i(k,q)
a996b8bdc0 Gael*0533              prof_interp_j(num_file,ProfNo_tile,q,bi,bj)
ba63501b4c Gael*0534      &       =tmp_j(k,q)
                0535              tmp_sum_weights=tmp_sum_weights+tmp_weights(k,q)
d5aa75d39a Jean*0536 c more test of the inputs: is the offline-computed
ba63501b4c Gael*0537 c interpolation information consistent (self and with grid)
                0538        if ( (tmp_i(k,q).LT.0).OR.(tmp_j(k,q).LT.0)
                0539      & .OR.(tmp_i(k,q).GT.sNx+1).OR.(tmp_j(k,q).GT.sNy+1) ) then
f0e4bffe35 Gael*0540           WRITE(msgBuf,'(4A)')
                0541      &     'PROFILES_INIT_FIXED: file ', profilesfile(1:IL),
                0542      &     '.nc includes inconsistent interpolation ',
                0543      &     'points (profilesDoGenGrid; out of tile)'
                0544           CALL PRINT_ERROR( msgBuf, myThid)
c9c84c2afb Gael*0545           stopGenericGrid=1
ba63501b4c Gael*0546        endif
4591803d6b Gael*0547 #ifdef ALLOW_PROFILES_EXCLUDE_CORNERS
c9c84c2afb Gael*0548        if ( tmp_weights(k,q) .NE. 0. _d 0) then
ba63501b4c Gael*0549        if ( ((tmp_i(k,q).EQ.0).AND.(tmp_j(k,q).EQ.0))
                0550      & .OR.((tmp_i(k,q).EQ.sNx+1).AND.(tmp_j(k,q).EQ.sNy+1))
                0551      & .OR.((tmp_i(k,q).EQ.0).AND.(tmp_j(k,q).EQ.sNy+1))
                0552      & .OR.((tmp_i(k,q).EQ.sNx+1).AND.(tmp_j(k,q).EQ.0)) ) then
f0e4bffe35 Gael*0553           WRITE(msgBuf,'(4A)')
                0554      &     'PROFILES_INIT_FIXED: file ', profilesfile(1:IL),
                0555      &     '.nc includes inconsistent interpolation ',
                0556      &     'points (profilesDoGenGrid; using overlap corners)'
                0557           CALL PRINT_ERROR( msgBuf, myThid)
c9c84c2afb Gael*0558           stopGenericGrid=1
ba63501b4c Gael*0559        endif
                0560        endif
4591803d6b Gael*0561 #endif /* ALLOW_PROFILES_EXCLUDE_CORNERS */
c9c84c2afb Gael*0562        if ( (tmp_weights(k,q).LT.0. _d 0).OR.
                0563      &    (tmp_weights(k,q).GT.1. _d 0) ) then
f0e4bffe35 Gael*0564           WRITE(msgBuf,'(4A)')
                0565      &     'PROFILES_INIT_FIXED: file ', profilesfile(1:IL),
                0566      &     '.nc includes inconsistent interpolation ',
                0567      &     'weights (profilesDoGenGrid; sum oustide 0-1)'
                0568           CALL PRINT_ERROR( msgBuf, myThid)
c9c84c2afb Gael*0569           stopGenericGrid=1
ba63501b4c Gael*0570        endif
                0571 
                0572        enddo
                0573 
c9c84c2afb Gael*0574        if ( abs(tmp_sum_weights -1. _d 0 ) .GT. 0.0001 _d 0) then
f0e4bffe35 Gael*0575           WRITE(msgBuf,'(4A)')
                0576      &     'PROFILES_INIT_FIXED: file ', profilesfile(1:IL),
                0577      &     '.nc includes inconsistent interpolation ',
                0578      &     'weights (profilesDoGenGrid; dont add up to 1)'
                0579           CALL PRINT_ERROR( msgBuf, myThid)
c9c84c2afb Gael*0580           stopGenericGrid=1
ba63501b4c Gael*0581        endif
                0582 
a996b8bdc0 Gael*0583          prof_ind_glob(num_file,ProfNo_tile,bi,bj)=k+1000*(kk-1)
ba63501b4c Gael*0584 
                0585        endif
                0586        endif
f0e4bffe35 Gael*0587        endif   !if (.NOT.profilesDoGenGrid) then
39ce977435 Gael*0588 
                0589 c ==============================================================================
                0590 
                0591 c check that maximum size was not reached:
a996b8bdc0 Gael*0592        if (ProfNo_tile.GE.NOBSGLOB) then
39ce977435 Gael*0593          WRITE(msgBuf,'(3A)')
                0594      &    'PROFILES_INIT_FIXED: file ', profilesfile(1:IL),
                0595      &    '.nc was not read properly (increase NOBSGLOB).'
                0596          CALL PRINT_ERROR( msgBuf, myThid)
                0597          stopProfiles=1
                0598        endif
                0599 
c9c84c2afb Gael*0600       endif    !if ( stopProfiles .EQ. 0) then
f0e4bffe35 Gael*0601       enddo    !do k=1,min(1000,ProfNo(num_file,bi,bj)-1000*(kk-1))
                0602       endif    !if (min(ProfNo(num_file,bi,bj), 1000...
                0603       enddo    !do kk=1,profno_div1000+1
d5aa75d39a Jean*0604 
a996b8bdc0 Gael*0605       ProfNo(num_file,bi,bj)=ProfNo_tile
                0606 
                0607       write(msgbuf,'(a,i9)')
                0608      &   '  # of profiles with erroneous HHMMSS values =',
                0609      &   ProfNo_hh
                0610       call print_message(
                0611      &   msgbuf, standardmessageunit, SQUEEZE_RIGHT , mythid)
ba63501b4c Gael*0612 
64278b1175 Jean*0613       write(msgbuf,'(a,i9)')
a996b8bdc0 Gael*0614      &   '  # of profiles within tile and time period  =',
                0615      &   ProfNo(num_file,bi,bj)
ba63501b4c Gael*0616       call print_message(
                0617      &   msgbuf, standardmessageunit, SQUEEZE_RIGHT , mythid)
                0618 
                0619 c6) available variablesin the data set
6a770e0a24 Patr*0620 
ea4d09597a Gael*0621       do k=1,NVARMAX
3c8dcfdea9 Gael*0622         prof_num_var_cur(num_file,k,bi,bj)=0
6a770e0a24 Patr*0623       enddo
71a5587721 Gael*0624       prof_num_var_tot(num_file,bi,bj)=0
6a770e0a24 Patr*0625 
3c8dcfdea9 Gael*0626       do k=1,NVARMAX
cf16ba6028 Gael*0627         JL  = ILNBLNK( prof_names(num_file,k) )
                0628         err = NF_INQ_VARID(fid,prof_names(num_file,k)(1:JL), varid1 )
3c8dcfdea9 Gael*0629         if (err.EQ.NF_NOERR) then
                0630           vec_quantities(num_file,k,bi,bj)=.TRUE.
                0631           prof_num_var_tot(num_file,bi,bj)=
                0632      &     prof_num_var_tot(num_file,bi,bj)+1
                0633           prof_num_var_cur(num_file,k,bi,bj)=
                0634      &     prof_num_var_tot(num_file,bi,bj)
                0635         else
                0636           vec_quantities(num_file,k,bi,bj)=.FALSE.
                0637         endif
                0638       enddo
6a770e0a24 Patr*0639 
38287224dd Gael*0640       do k=1,NVARMAX
                0641         if (vec_quantities(num_file,k,bi,bj)) then
                0642           KL  = ILNBLNK( prof_names(num_file,k) )
                0643           JL  = ILNBLNK( prof_namesmod(num_file,k) )
                0644           if (prof_namesmod(num_file,k).EQ.'pTracer') then
                0645       write(msgbuf,'(a,I3,5a,I3)') '  variable #',k,' is ' ,
                0646      & prof_names(num_file,k)(1:KL),' and ',
                0647      & prof_namesmod(num_file,k)(1:JL),' #',
                0648      & prof_itracer(num_file,k)
                0649           else
a996b8bdc0 Gael*0650       write(msgbuf,'(a,I3,4a)') '  variable #',k,
                0651      & ' is            ' ,
38287224dd Gael*0652      & prof_names(num_file,k)(1:KL),' and ',
                0653      & prof_namesmod(num_file,k)(1:JL)
                0654           endif
                0655           call print_message(msgbuf,
                0656      &       standardmessageunit, SQUEEZE_RIGHT , mythid)
                0657         endif
                0658       enddo
6a770e0a24 Patr*0659 
                0660 C===========================================================
                0661 c create files for model counterparts to observations
                0662 C===========================================================
                0663 
d5aa75d39a Jean*0664            if (ProfNo(num_file,bi,bj).GT.0) then
6a770e0a24 Patr*0665          iG=bi+(myXGlobalLo-1)/sNx ! Kludge until unstructered tiles
                0666          jG=bj+(myYGlobalLo-1)/sNy ! Kludge until unstructered tiles
                0667 
1ff0163ead Gael*0668       JL  = ILNBLNK( profilesDir )
                0669 
f0e4bffe35 Gael*0670       if (profilesDoNcOutput) then
d5aa75d39a Jean*0671 
de57a2ec4b Mart*0672       write(fnameequinc,'(3a,i3.3,a,i3.3,a)')
1ff0163ead Gael*0673      & profilesDir(1:JL),profilesfile(1:IL),'.',iG,'.',jG,'.equi.nc'
de57a2ec4b Mart*0674       write(adfnameequinc,'(4a,i3.3,a,i3.3,a)')
1ff0163ead Gael*0675      & profilesDir(1:JL),'ad',
6a770e0a24 Patr*0676      & profilesfile(1:IL),'.',iG,'.',jG,'.equi.nc'
                0677 
                0678       inquire( file=fnameequinc, exist=exst )
                0679       if (.NOT.exst) then
b2a948f981 Gael*0680         call profiles_init_ncfile(num_file,
                0681      &   fiddata(num_file,bi,bj),fnameequinc,
                0682      &   fidforward(num_file,bi,bj),ProfNo(num_file,bi,bj),
                0683      &   ProfDepthNo(num_file,bi,bj),
                0684      &   bi,bj,myThid)
6a770e0a24 Patr*0685       else
b2a948f981 Gael*0686         err = NF_OPEN(fnameequinc,NF_WRITE,fidforward(num_file,bi,bj))
6a770e0a24 Patr*0687       endif
b2a948f981 Gael*0688 #ifdef ALLOW_ADJOINT_RUN
                0689       inquire( file=adfnameequinc, exist=exst )
                0690       if (.NOT.exst) then
                0691         call profiles_init_ncfile(num_file,fiddata(num_file,bi,bj),
                0692      &   adfnameequinc, fidadjoint(num_file,bi,bj),
                0693      &   ProfNo(num_file,bi,bj),ProfDepthNo(num_file,bi,bj),
                0694      & bi,bj, myThid)
                0695       else
                0696         err = NF_OPEN(adfnameequinc,NF_WRITE,fidadjoint(num_file,bi,bj))
                0697       endif
                0698 #endif
6a770e0a24 Patr*0699       else
                0700 
de57a2ec4b Mart*0701       write(fnameequinc,'(3a,i3.3,a,i3.3,a)')
1ff0163ead Gael*0702      & profilesDir(1:JL),profilesfile(1:IL),'.',iG,'.',jG,'.equi.data'
de57a2ec4b Mart*0703       write(adfnameequinc,'(4a,i3.3,a,i3.3,a)')
1ff0163ead Gael*0704      & profilesDir(1:JL),'ad',
64e4a6baa3 Jean*0705      & profilesfile(1:IL),'.',iG,'.',jG,'.equi.data'
6a770e0a24 Patr*0706 
                0707       inquire( file=fnameequinc, exist=exst )
4caf2dc194 Gael*0708 #ifdef PROFILES_USE_MDSFINDUNITS
                0709       call MDSFINDUNIT( fidforward(num_file,bi,bj) , mythid )
                0710 #else
                0711       call PROFILES_FINDUNIT( fidforward(num_file,bi,bj) , mythid )
                0712 #endif
6a770e0a24 Patr*0713       if (.NOT.exst) then
b2a948f981 Gael*0714         call profiles_init_ncfile(num_file,fiddata(num_file,bi,bj),
                0715      &   fnameequinc,fidforward(num_file,bi,bj),
                0716      &   ProfNo(num_file,bi,bj),ProfDepthNo(num_file,bi,bj),
                0717      &   bi,bj,myThid)
6a770e0a24 Patr*0718       else
b2a948f981 Gael*0719          open( fidforward(num_file,bi,bj),file=fnameequinc,
                0720      &   form ='unformatted',status='unknown', access='direct',
                0721      &   recl=  (ProfDepthNo(num_file,bi,bj)+1)*WORDLENGTH*2 )
                0722       endif
51c7f9a83e Gael*0723 #ifdef ALLOW_ADJOINT_RUN
b2a948f981 Gael*0724       inquire( file=adfnameequinc, exist=exst )
4caf2dc194 Gael*0725 #ifdef PROFILES_USE_MDSFINDUNITS
                0726       call MDSFINDUNIT( fidadjoint(num_file,bi,bj) , mythid )
                0727 #else
                0728       call PROFILES_FINDUNIT( fidadjoint(num_file,bi,bj) , mythid )
                0729 #endif
b2a948f981 Gael*0730       if (.NOT.exst) then
                0731         call profiles_init_ncfile(num_file,fiddata(num_file,bi,bj),
                0732      &   adfnameequinc, fidadjoint(num_file,bi,bj),
                0733      &   ProfNo(num_file,bi,bj),ProfDepthNo(num_file,bi,bj),
                0734      &   bi,bj, myThid)
                0735       else
                0736          open( fidadjoint(num_file,bi,bj),file=adfnameequinc,
                0737      &   form ='unformatted',status='unknown', access='direct',
                0738      &   recl=  (ProfDepthNo(num_file,bi,bj)+1)*WORDLENGTH*2 )
6a770e0a24 Patr*0739       endif
b2a948f981 Gael*0740 #endif
6a770e0a24 Patr*0741 
                0742       endif
                0743 
                0744            endif
                0745 
                0746 C===========================================================
                0747       else
71a5587721 Gael*0748       ProfNo(num_file,bi,bj)=0
ea4d09597a Gael*0749       do k=1,NVARMAX
71a5587721 Gael*0750       prof_num_var_cur(num_file,k,bi,bj)=0
                0751       vec_quantities(num_file,k,bi,bj)=.FALSE.
6a770e0a24 Patr*0752       enddo
71a5587721 Gael*0753       prof_num_var_tot(num_file,bi,bj)=0
6a770e0a24 Patr*0754       do k=1,NOBSGLOB
c9c84c2afb Gael*0755       prof_time(num_file,k,bi,bj)=-999. _d 0
                0756       prof_lon(num_file,k,bi,bj)=-999. _d 0
                0757       prof_lat(num_file,k,bi,bj)=-999. _d 0
                0758       prof_ind_glob(num_file,k,bi,bj)=0
6b2230d510 Ou W*0759 #ifdef ALLOW_PROFILES_SAMPLESPLIT_COST
                0760       prof_ind_avgbin(num_file,k,bi,bj)=-999
                0761 #endif
ba63501b4c Gael*0762       do q = 1,NUM_INTERP_POINTS
38287224dd Gael*0763          prof_interp_i(num_file,k,q,bi,bj) = 1
                0764          prof_interp_j(num_file,k,q,bi,bj) = 1
                0765          prof_interp_weights(num_file,k,q,bi,bj) = 0. _d 0
ba63501b4c Gael*0766       enddo
c9c84c2afb Gael*0767       prof_interp_xC11(num_file,k,bi,bj)=-999. _d 0
                0768       prof_interp_yC11(num_file,k,bi,bj)=-999. _d 0
                0769       prof_interp_xCNINJ(num_file,k,bi,bj)=-999. _d 0
                0770       prof_interp_yCNINJ(num_file,k,bi,bj)=-999. _d 0
6a770e0a24 Patr*0771       enddo
                0772 
                0773       endif !if (IL.NE.0) then
                0774       enddo !      do num_file=1,NFILESPROFMAX
ba63501b4c Gael*0775 
6b2230d510 Ou W*0776 #ifdef ALLOW_PROFILES_SAMPLESPLIT_COST
                0777 C Find the unique depth levels from all profile datasets   
                0778 C initialize prof_depth_comb
                0779       if(bi.EQ.1.and.bj.EQ.1)then
                0780          NLEVELCOMB = 0
                0781          NLEVELCOMBRL = NLEVELCOMB
                0782       endif
                0783       do m=1,NLEVELCOMBMAX
                0784          prof_depth_comb(m,bi,bj)=-999. _d 0
                0785       enddo
                0786 
                0787       m = 1
                0788       do num_file=1,NFILESPROFMAX
                0789        do K=1,ProfDepthNo(num_file,bi,bj)
                0790 
                0791           if(m.EQ.1) then
                0792            prof_depth_comb(m,bi,bj) = prof_depth(num_file, k,bi,bj)
                0793            m = m + 1
                0794           else
                0795 C sort
                0796            do l=1,NLEVELCOMBMAX-1
                0797             if(prof_depth_comb(l,bi,bj) .ne. -999. _d 0) then
                0798 
                0799               if(prof_depth(num_file, k,bi,bj).lt.
                0800      &           prof_depth_comb(l,bi,bj).and.
                0801      &           l.EQ.1)  then
                0802                  prof_depth_comb(NLEVELCOMBMAX,bi,bj) =
                0803      &            prof_depth_comb(l,bi,bj)
                0804                  prof_depth_comb(l,bi,bj)=
                0805      &            prof_depth(num_file, k,bi,bj)
                0806                  do il = NLEVELCOMBMAX-1, l+2,-1
                0807                     prof_depth_comb(il,bi,bj)=
                0808      &            prof_depth_comb(il-1,bi,bj)
                0809                  enddo
                0810                  prof_depth_comb(l+1,bi,bj)=
                0811      &            prof_depth_comb(NLEVELCOMBMAX,bi,bj)
                0812               else if(prof_depth(num_file, k,bi,bj).gt.
                0813      &           prof_depth_comb(l,bi,bj).and.
                0814      &           prof_depth(num_file, k,bi,bj).lt.
                0815      &           prof_depth_comb(l+1,bi,bj))  then
                0816 
                0817                  prof_depth_comb(NLEVELCOMBMAX,bi,bj) =
                0818      &            prof_depth_comb(l+1,bi,bj)
                0819                  prof_depth_comb(l+1,bi,bj)=
                0820      &            prof_depth(num_file, k,bi,bj)
                0821                  do il = NLEVELCOMBMAX-1, l+3,-1
                0822                     prof_depth_comb(il,bi,bj)=
                0823      &            prof_depth_comb(il-1,bi,bj)
                0824                  enddo
                0825                  prof_depth_comb(l+2,bi,bj)=
                0826      &            prof_depth_comb(NLEVELCOMBMAX,bi,bj)
                0827               else if ( prof_depth(num_file, k,bi,bj).gt.
                0828      &           prof_depth_comb(l,bi,bj).and.
                0829      &           prof_depth_comb(l+1,bi,bj).eq.-999. _d 0)  then
                0830                  prof_depth_comb(l+1,bi,bj) =
                0831      &              prof_depth(num_file, k,bi,bj)
                0832               endif
                0833              endif
                0834            enddo
                0835 
                0836           endif
                0837           if(m.GE.NLEVELCOMBMAX-2)then
                0838             WRITE(msgBuf,'(A)')
                0839      &      'increase NLEVELCOMBMAX'
                0840             CALL PRINT_ERROR( msgBuf, myThid)
                0841            endif
                0842        enddo ! do K=1,ProfDepthNo(num_file,bi,bj)
                0843       enddo ! do num_file=1,NFILESPROFMAX
                0844       prof_depth_comb(NLEVELCOMBMAX,bi,bj) = -999. _d 0
                0845 
                0846 C diagnostics output
                0847       do m=1,NLEVELCOMBMAX
                0848          if(prof_depth_comb(m,bi,bj) .GE. 0. _d 0
                0849      &     .AND. NLEVELCOMB.LT.m)then
                0850            NLEVELCOMB = m
                0851            if(m.GE.NLEVELCOMBMAX-2)then
                0852             WRITE(msgBuf,'(A,2i6)')
                0853      &      'increase NLEVELCOMBMAX: m,NLEVELCOMBMA  ',
                0854      &      m, NLEVELCOMBMAX
                0855             CALL PRINT_ERROR( msgBuf, myThid)
                0856            endif
                0857          endif
                0858       enddo
                0859       WRITE(msgBuf,'(A, i6,d20.5)')
                0860      &      'NLEVELCOMB = ', NLEVELCOMB
                0861       call print_message( msgbuf, standardmessageunit,
                0862      &                    SQUEEZE_RIGHT , mythid)
                0863 #endif
                0864 
6a770e0a24 Patr*0865 C===========================================================
ba63501b4c Gael*0866 C error cases:
                0867 C===========================================================
                0868 
                0869 c1) you want to provide interpolation information
                0870 
c9c84c2afb Gael*0871        if ( stopGenericGrid.EQ.2) then
ba63501b4c Gael*0872          iG=bi+(myXGlobalLo-1)/sNx ! Kludge until unstructered tiles
                0873          jG=bj+(myYGlobalLo-1)/sNy ! Kludge until unstructered tiles
                0874 cgf XC grid
                0875        call MDSFINDUNIT( fid , mythid )
de57a2ec4b Mart*0876        write(fnameequinc,'(a,i3.3,a,i3.3,a,i4.4,a,i4.4,a)')
ba63501b4c Gael*0877      & 'profilesXCincl1PointOverlap.',iG,'.',jG,'.',sNx,'.',sNy,'.data'
                0878          k=MDS_RECLEN(64,(sNx+2)*(sNy+2),mythid)
                0879             WRITE(standardMessageUnit,'(A,/,2A)')
                0880      & 'PROFILES_INIT_FIXED: creating grid from profiles; file:',
                0881      & fnameequinc
d5aa75d39a Jean*0882        open( fid, file= fnameequinc, form ='unformatted',
ba63501b4c Gael*0883      &      status='unknown',access='direct', recl= k)
                0884         DO m=0,sNy+1
                0885          DO l=0,sNx+1
                0886         xy_buffer_r8(l,m)=xC(l,m,bi,bj)
                0887          ENDDO
                0888         ENDDO
                0889 #ifdef _BYTESWAPIO
                0890             call MDS_BYTESWAPR8((sNx+2)*(sNy+2),xy_buffer_r8)
                0891 #endif
                0892        write(fid,rec=1) xy_buffer_r8
                0893        close(fid)
                0894 cgf YC grid
                0895        call MDSFINDUNIT( fid , mythid )
de57a2ec4b Mart*0896        write(fnameequinc,'(a,i3.3,a,i3.3,a,i4.4,a,i4.4,a)')
ba63501b4c Gael*0897      & 'profilesYCincl1PointOverlap.',iG,'.',jG,'.',sNx,'.',sNy,'.data'
                0898          k=MDS_RECLEN(64,(sNx+2)*(sNy+2),mythid)
                0899             WRITE(standardMessageUnit,'(A,/,A)')
                0900      & 'PROFILES_INIT_FIXED: creating grid from profiles; file:',
                0901      & fnameequinc
d5aa75d39a Jean*0902        open( fid, file= fnameequinc, form ='unformatted',
ba63501b4c Gael*0903      & status='unknown', access='direct', recl= k)
                0904         DO m=0,sNy+1
                0905          DO l=0,sNx+1
                0906                 xy_buffer_r8(l,m)=yC(l,m,bi,bj)
                0907          ENDDO
                0908         ENDDO
                0909 #ifdef _BYTESWAPIO
                0910             call MDS_BYTESWAPR8((sNx+2)*(sNy+2),xy_buffer_r8)
                0911 #endif
                0912        write(fid,rec=1) xy_buffer_r8
                0913        close(fid)
f0e4bffe35 Gael*0914 
                0915        WRITE(msgBuf,'(3A)')
                0916      & 'PROFILES_INIT_FIXED : ',
39ce977435 Gael*0917      & 'when using profilesDoGenGrid ',
f0e4bffe35 Gael*0918      & 'you have to provide interpolation coeffs etc. '
                0919        CALL PRINT_ERROR( msgBuf, myThid)
                0920        WRITE(msgBuf,'(2A)')
                0921      & 'and some of your nc files dont have them. ',
                0922      & 'You could use profiles_prep_mygrid.m and/or'
                0923        CALL PRINT_ERROR( msgBuf, myThid)
                0924        WRITE(msgBuf,'(A)')
                0925      & 'use the grid info in profiles*incl1PointOverlap*data'
                0926        CALL PRINT_ERROR( msgBuf, myThid)
c9c84c2afb Gael*0927        stopProfiles=1
ba63501b4c Gael*0928 
                0929       endif
                0930 
71a5587721 Gael*0931       ENDDO
                0932       ENDDO
                0933 
6b2230d510 Ou W*0934 #ifdef ALLOW_PROFILES_SAMPLESPLIT_COST
                0935       NLEVELCOMBRL = NLEVELCOMB
                0936       _GLOBAL_MAX_RL( NLEVELCOMBRL, myThid )
                0937       NLEVELCOMB = NLEVELCOMBRL
                0938 #endif
                0939 
ba63501b4c Gael*0940       _END_MASTER( mythid )
                0941       _BARRIER
                0942 
                0943 c2) stop after other kind of errors
c9c84c2afb Gael*0944       CALL GLOBAL_SUM_INT( stopProfiles , myThid )
                0945       if ( stopProfiles.GE.1) then
b00d6c1700 Gael*0946         CALL ALL_PROC_DIE( myThid )
                0947         STOP 'ABNORMAL END: S/R PROFILES_INIT_FIXED'
ba63501b4c Gael*0948       endif
39ce977435 Gael*0949 
c9c84c2afb Gael*0950       CALL GLOBAL_SUM_INT( stopGenericGrid , myThid )
                0951       if ( stopGenericGrid.GE.1) then
b00d6c1700 Gael*0952         CALL ALL_PROC_DIE( myThid )
                0953         STOP 'ABNORMAL END: S/R PROFILES_INIT_FIXED'
ba63501b4c Gael*0954       endif
                0955 
f0e4bffe35 Gael*0956       write(msgbuf,'(a)') ' '
                0957       call print_message( msgbuf, standardmessageunit,
                0958      &                    SQUEEZE_RIGHT , mythid)
                0959       write(msgbuf,'(a)')
                0960      &'// ======================================================='
                0961       call print_message( msgbuf, standardmessageunit,
                0962      &                    SQUEEZE_RIGHT , mythid)
                0963       write(msgbuf,'(a)')
                0964      &'// insitu profiles model sampling >>> END <<<'
                0965       call print_message( msgbuf, standardmessageunit,
                0966      &                    SQUEEZE_RIGHT , mythid)
                0967       write(msgbuf,'(a)')
                0968      &'// ======================================================='
                0969       call print_message( msgbuf, standardmessageunit,
                0970      &                    SQUEEZE_RIGHT , mythid)
                0971       write(msgbuf,'(a)') ' '
                0972       call print_message( msgbuf, standardmessageunit,
                0973      &                    SQUEEZE_RIGHT , mythid)
                0974 
6e4c90fea3 Patr*0975 #endif
6a770e0a24 Patr*0976 
d5aa75d39a Jean*0977       RETURN
6a770e0a24 Patr*0978       END