Index: /issm/trunk-jpl/test/NightlyRun/test4002.m
===================================================================
--- /issm/trunk-jpl/test/NightlyRun/test4002.m	(revision 22671)
+++ /issm/trunk-jpl/test/NightlyRun/test4002.m	(revision 22672)
@@ -573,7 +573,4 @@
         !\cp ../MITgcm/build/mitgcmuv run
         !\cp ../MITgcm/input/* run
-
-        %load data: 
-        loaddata(org,'Parameters');
 
         % initial salinity
@@ -845,9 +842,146 @@
 	loaddata(org,'Parameters');
 	loaddata(org,'ParameterizeIce');
-
-	if exist ('run')
-		!\rm -rf run/test4002*
-	end
-	
+	loaddata(org,'Bathymetry');
+	loaddata(org,'IceSheetGeometry');
+        endtime = round(MITgcmDeltaT * ...
+         floor(time_step*y2s*async_step_MITgcm_multiplier/MITgcmDeltaT));
+
+        % {{{ prepare MITgcm 
+        % rename previous run directory and create new one
+        if exist ('run.old')
+            !\rm -rf run.old
+        end
+        if exist ('run')
+            !\mv run run.old
+        end
+        !\mkdir run
+        !\cp ../MITgcm/build/mitgcmuv run
+        !\cp ../MITgcm/input/* run
+
+        % initial salinity
+        S=iniSalt*ones(Nx,Ny,Nz);
+        writebin('run/Salt.bin',S);
+
+        % initial temperature
+        T=iniTheta*ones(Nx,Ny,Nz);
+        writebin('run/Theta.bin',T);
+
+        % initial velocity
+        Z=zeros(Nx,Ny,Nz);
+        writebin('run/Uvel.bin',Z);
+        writebin('run/Vvel.bin',Z);
+
+        % initial sea surface height
+        Z=zeros(Nx,Ny);
+        writebin('run/Etan.bin',Z);
+
+        % salinity boundary conditions
+        S=obcSalt*ones(Ny,Nz);
+        thk=delZ*ones(Nz,1);
+        bot=cumsum(thk);
+        ik=find(bot<=mlDepth);
+        S(:,ik)=mlSalt;
+        writebin('run/OBs.bin',S);
+
+        % temperature boundary conditions
+        T=obcTheta*ones(Ny,Nz);
+        T(:,ik)=mlTheta;
+        writebin('run/OBt.bin',T);
+
+        % zonal velocity boundary conditions
+        U=obcUvel*ones(Ny,Nz);
+        writebin('run/OBu.bin',U);
+
+        % zero boundary conditions
+        Z=zeros(Ny,Nz);
+        writebin('run/zeros.bin',Z);
+
+        % build parameter file data.obcs
+        fidi=fopen('../MITgcm/input/data.obcs','r');
+        fido=fopen('run/data.obcs','w');
+        tline = fgetl(fidi);
+        fprintf(fido,'%s\n',tline);
+        while 1
+            tline = fgetl(fidi);
+            if ~ischar(tline), break, end
+            %do the change here: 
+            if strcmpi(tline,' OB_Iwest = 40*1,'),
+                fprintf(fido,'%s%i%s\n',' OB_Iwest = ',Ny,'*1,');
+                continue;
+            end
+            if strcmpi(tline,' OB_Ieast = 40*-1,'),
+                fprintf(fido,'%s%i%s\n',' OB_Ieast = ',Ny,'*-1,');
+                continue;
+            end
+            fprintf(fido,'%s\n',tline);
+        end
+        %close  files
+        fclose(fidi);
+        fclose(fido);
+
+        %save bathymetry and bedrock in run directory
+        writebin('run/bathymetry.bin',bathymetry);
+        writebin('run/icetopo.bin',draft);
+        % }}}
+        % {{{ generate MITgcm parameter file data 
+        fidi=fopen('../MITgcm/input/data','r');
+        fido=fopen('run/data','w');
+        tline = fgetl(fidi);
+        fprintf(fido,'%s\n',tline);
+        while 1
+            tline = fgetl(fidi);
+            if ~ischar(tline), break, end
+            %do the change here: 
+            if strcmpi(tline,' xgOrigin = 0.0,'),
+                fprintf(fido,'%s%i%s\n',' xgOrigin = ',xgOrigin,',');
+                continue;
+            end
+            if strcmpi(tline,' ygOrigin = -80.0,'),
+                fprintf(fido,'%s%i%s\n',' ygOrigin = ',ygOrigin,',');
+                continue;
+            end
+            if strcmpi(tline,' delX = 20*0.25,'),
+                fprintf(fido,'%s%i*%g%s\n',' delX = ',Nx,dLong,',');
+                continue;
+            end
+            if strcmpi(tline,' delY = 20*0.25,'),
+                fprintf(fido,'%s%i*%g%s\n',' delY = ',Ny,dLat,',');
+                continue;
+            end
+            if strcmpi(tline,' delZ = 30*30.0,'),
+                fprintf(fido,'%s%i*%g%s\n',' delZ = ',Nz,delZ,',');
+                continue;
+            end
+            if strcmpi(tline,' endTime=2592000.,'),
+                fprintf(fido,'%s%i%s\n',' endTime= ',endtime,',');
+                continue;
+            end
+            if strcmpi(tline,' deltaT=1200.0,'),
+                fprintf(fido,'%s%i%s\n',' deltaT= ',MITgcmDeltaT,',');
+                continue;
+            end
+            if strcmpi(tline,' pChkptFreq=2592000.,'),
+                fprintf(fido,'%s%i%s\n',' pChkptFreq= ',endtime,',');
+                continue;
+            end
+            if strcmpi(tline,' taveFreq=2592000.,'),
+                fprintf(fido,'%s%i%s\n',' taveFreq= ',endtime,',');
+                continue;
+            end
+            if strcmpi(tline,' rhoConst=1030.,'),
+                fprintf(fido,'%s%i%s\n',' rhoConst= ',rho_water,',');
+                continue;
+            end
+            if strcmpi(tline,' rhoNil=1030.,'),
+                fprintf(fido,'%s%i%s\n',' rhoNil= ',rho_water,',');
+                continue;
+            end
+            fprintf(fido,'%s\n',tline);
+        end
+        %close  files
+        fclose(fidi);
+        fclose(fido);
+        % }}}
+	        
 	md.transient.isoceancoupling=1;
 	md.transient.isgroundingline=1;
