#include "CPP_OPTIONS.h" subroutine smooth_ensinv (mythid) IMPLICIT NONE #include "SIZE.h" #include "EEPARAMS.h" #include "GRID.h" #include "PARAMS.h" c#include "tamc.h" #include "smooth.h" #ifdef ALLOW_PROFILES #include "profiles.h" #include "netcdf.inc" #endif #include "DYNVARS.h" c#include "adcommon.h" #include "optim.h" #include "ecco_cost.h" integer myThid #ifdef ALLOW_SMOOTH_ENSINV double precision adsalt(1-olx:snx+olx,1-oly:sny+oly,nr,nsx,nsy) double precision adtheta(1-olx:snx+olx,1-oly:sny+oly,nr,nsx,nsy) common /addynvars_r/ adtheta, adsalt character*( 80) fnamegeneric integer i,j,k,bi,bj integer itlo,ithi integer jtlo,jthi _RL port_rand, port_rand_norm integer iloop,iii,ii,nbRand integer nbRandPrev integer choice_bckgd _RL mytime integer myiter integer num_file,num_var,prof_num integer iG,jG,err,fid _RL tmp_lon _RL prof_traj1D(NLEVELMAX), prof_traj1D_mean(NLEVELMAX) _RL prof_data1D(NLEVELMAX), prof_weights1D(NLEVELMAX) integer nbt_in jtlo = mybylo(mythid) jthi = mybyhi(mythid) itlo = mybxlo(mythid) ithi = mybxhi(mythid) write(fnamegeneric(1:80),'(1a)') & 'wvar_T_4_50levels.bin' call mdsreadfield(fnamegeneric,32,'RL',nR, & wthetalev,1,mythid) _EXCH_XYZ_RL( wthetalev, mythid ) write(fnamegeneric(1:80),'(1a)') & 'wvar_S_4_50levels.bin' call mdsreadfield(fnamegeneric,32,'RL',nR, & wsaltlev,1,mythid) _EXCH_XYZ_RL( wsaltlev, mythid ) C initialize the random number generator: theta(1,1,1,1,1)=port_rand(1) nbRand=200 nbRandPrev=119 choice_bckgd=0 C load the background fields: if (choice_bckgd.GT.0) then write(fnamegeneric(1:80),'(1a)') & 'wc01_3Densinv.bckgd' call mdsreadfield(fnamegeneric,32,'RL',nR,wtheta2, & 1,mythid) call mdsreadfield(fnamegeneric,32,'RL',nR,wsalt2, & 2,mythid) endif C main loop: DO ii=1,nbRand print*,"ii",ii C define a perturbation: DO bj=jtlo,jthi DO bi=itlo,ithi DO k=1,Nr DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx theta(i,j,k,bi,bj)=0. salt(i,j,k,bi,bj)=0. if (maskC(i,j,k,bi,bj).NE.0) then theta(i,j,k,bi,bj)=port_rand_norm() salt(i,j,k,bi,bj)=port_rand_norm() if (choice_bckgd.EQ.1) then theta(i,j,k,bi,bj)=theta(i,j,k,bi,bj)+wtheta2(i,j,k,bi,bj) salt(i,j,k,bi,bj)=salt(i,j,k,bi,bj)+wsalt2(i,j,k,bi,bj) endif endif ENDDO ENDDO ENDDO ENDDO ENDDO _EXCH_XYZ_RL ( theta, myThid ) _EXCH_XYZ_RL ( salt, myThid ) if (ii.GT.nbRandPrev) then write(fnamegeneric(1:80),'(1a,i3.3)') & 'wc01_3Densinv.',ii C the smoothing itself: call smooth_correl3D(theta,1,mythid) call smooth_correl3D(salt,1,mythid) C scale the variance DO bj=jtlo,jthi DO bi=itlo,ithi DO k=1,Nr DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx if (maskC(i,j,k,bi,bj).NE.0) then theta(i,j,k,bi,bj)=theta(i,j,k,bi,bj) & * wthetaLev(i,j,k,bi,bj) * maskc(i,j,k,bi,bj) salt(i,j,k,bi,bj)=salt(i,j,k,bi,bj) & * wsaltLev(i,j,k,bi,bj) * maskc(i,j,k,bi,bj) if (choice_bckgd.EQ.2) then theta(i,j,k,bi,bj)=theta(i,j,k,bi,bj)+wtheta2(i,j,k,bi,bj) salt(i,j,k,bi,bj)=salt(i,j,k,bi,bj)+wsalt2(i,j,k,bi,bj) endif endif ENDDO ENDDO ENDDO ENDDO ENDDO _EXCH_XYZ_RL ( theta, myThid ) _EXCH_XYZ_RL ( salt, myThid ) c do the profiles_inloop operations do iloop=1,nTimeSteps mytime = startTime + float(iloop-1)*deltaTclock call profiles_inloop( mytime, mythid ) enddo c compute the cost and write to ad*equi.bin DO bj=jtlo,jthi DO bi=itlo,ithi do num_file=1,NFILESPROFMAX fid=fiddata(num_file,bi,bj) do prof_num=1,NOBSGLOB if (prof_num.LE.ProfNo(num_file,bi,bj)) then do num_var=1,NVARMAX do K=1,NLEVELMAX prof_traj1D(k)=0. prof_mask1D_cur(k,bi,bj)=0. prof_data1D(k)=0. prof_weights1D(k)=0. enddo if (vec_quantities(num_file,num_var,bi,bj).EQV..TRUE.) then do K=1,ProfDepthNo(num_file,bi,bj) prof_traj1D_mean(K)=0. enddo call active_read_profile(num_file, & ProfDepthNo(num_file,bi,bj),prof_traj1D,num_var, & prof_num,.false.,optimcycle,bi,bj,mythid, & profiles_dummy(num_file,num_var,bi,bj)) call profiles_readvector(num_file,num_var, & prof_ind_glob(num_file,prof_num,bi,bj), & ProfDepthNo(num_file,bi,bj),prof_data1D,bi,bj,myThid) call profiles_readvector(num_file,-num_var, & prof_ind_glob(num_file,prof_num,bi,bj), & ProfDepthNo(num_file,bi,bj), & prof_weights1D,bi,bj,myThid) do K=1,ProfDepthNo(num_file,bi,bj) prof_traj1D(K)= & prof_weights1D(K)*prof_mask1D_cur(K,bi,bj) & *(prof_traj1D(K)-prof_data1D(K)-prof_traj1D_mean(K)) & *(prof_traj1D(K)-prof_data1D(K)-prof_traj1D_mean(K)) enddo call adactive_read_profile( num_file,profdepthno(num_file,bi,bj), &num_var,prof_num, .false. ,optimcycle,bi,bj,mythid,prof_traj1d ) endif enddo !do num_var... endif !if (prof_num.LE.ProfNo(num_file,bi,bj)) then enddo !do prof_num=.. enddo !do num_file... ENDDO ENDDO c initialise adtheta/adsalt and do adinloop DO bj=jtlo,jthi DO bi=itlo,ithi DO k=1,Nr DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx adtheta(i,j,k,bi,bj)=0. adsalt(i,j,k,bi,bj)=0. ENDDO ENDDO ENDDO ENDDO ENDDO _EXCH_XYZ_RL ( adtheta, myThid ) _EXCH_XYZ_RL ( adsalt, myThid ) do iloop = ntimesteps, 1, -1 mytime = starttime+float(iloop-1)*deltatclock call adprofiles_inloop( mytime,mythid ) end do c co adexch and exch call adexch_xyz_rl( mythid,adtheta ) call adexch_xyz_rl( mythid,adsalt ) call exch_xyz_rl( adtheta,mythid ) call exch_xyz_rl( adsalt,mythid ) c do smooth_diff c nbt_in=wc01_nbt(smoothOpNbCur) c call smooth_diff3D(adtheta,nbt_in,mythid) c call smooth_diff3D(adsalt,nbt_in,mythid) call mdswritefield(fnamegeneric,32,.false.,'RL', & nR,theta,1,1,mythid) call mdswritefield(fnamegeneric,32,.false.,'RL', & nR,adtheta,2,1,mythid) call mdswritefield(fnamegeneric,32,.false.,'RL', & nR,salt,3,1,mythid) call mdswritefield(fnamegeneric,32,.false.,'RL', & nR,adsalt,4,1,mythid) endif ENDDO stop #endif / * ALLOW_SMOOTH_ENSINV * / end