FFT UPDATES

brayantz · Aug 30, 2012 · 14b717f · 14b717f
1 parent 8c572b5
commit 14b717f
Show file tree

Hide file tree

Showing 5 changed files with 41 additions and 24 deletions.
diff --git a/fft_ma_1d.m b/fft_ma_1d.m
@@ -56,7 +56,17 @@
 
 nx_c=nx*options.fac_x;
 
-% COVARIANCE MODEL
+
+
+%% REMOVE OLD COVARIANCE OF options.constant_C=0
+if (isfield(options,'constant_C'));
+    if options.constant_C==0;
+        try;options=rmfield(options,'C');end
+        try;options=rmfield(options,'fftC');end
+    end
+end
+
+%% COVARIANCE MODEL
 if (~isfield(options,'C'))&(~isfield(options,'fftC'));
     options.C=zeros(1,nx_c);
 

diff --git a/fft_ma_2d.m b/fft_ma_2d.m
@@ -126,7 +126,14 @@
 nx_c=nx+options.pad_x;
 x_all=[0:1:(nx_c-1)].*dx+x(1);
 y_all=[0:1:(ny_c-1)].*dy+y(1);
-
+
+%% REMOVE OLD COVARIANCE OF options.constant_C=0
+if (isfield(options,'constant_C'));
+    if options.constant_C==0;
+        try;options=rmfield(options,'C');end
+        try;options=rmfield(options,'fftC');end
+    end
+end
 
 %% SETUP  COVARIANCE MODEL
 if (~isfield(options,'C'))&(~isfield(options,'fftC'));

diff --git a/fft_ma_3d.m b/fft_ma_3d.m
@@ -157,6 +157,14 @@
 y_all=[0:1:(ny_c-1)].*dy;
 z_all=[0:1:(nz_c-1)].*dz;
 
+%% REMOVE OLD COVARIANCE OF options.constant_C=0
+if (isfield(options,'constant_C'));
+    if options.constant_C==0;
+        try;options=rmfield(options,'C');end
+        try;options=rmfield(options,'fftC');end
+    end
+end
+
 %% SETUP  COVARIANCE MODEL
 if (~isfield(options,'C'))&(~isfield(options,'fftC'));
 

diff --git a/krig_covar_lik.m b/krig_covar_lik.m
@@ -37,33 +37,22 @@
     Cm=precal_cov(pos_known,pos_known,V,options);
     %Cm=precal_cov(pos_known,pos_known,V);
 end
-% TMH 05/05/2009 :
-% CM IS NOT RIGHT IN CASE OF A NUGGET
-% Cm=precal_cov(pos_known,pos_known,V)
-% NUGGET IS NOT PRESENT
-
 
 dm=val_known(:,1)-m0;
 
-
-
-
 if method==1
     mgstat_verbose(sprintf('%s : Pardo-Iguzquiza likelihood',mfilename),1);
     % FROM Pardo-Iguzquiza, 1998, Math Geol, 30(1), EQN 8.
     % Kitanidis (1983), Water Res. Res. 19(4), EQN 29.
     sample_size=nknown;
     sill=sum([V.par1]);
 
-    %Cd=eye(size(Cm))+0.00001*sill; % FOR BETTER PERFORMANCE
-    %Cd=eye(size(Cm)).*0.00001*sill; % FOR BETTER PERFORMANCE
     Cd=eye(size(Cm)).*0.001*sill; % FOR BETTER PERFORMANCE
     Cm=Cm+Cd;
 
     Q=Cm./sill;
     d_val=val_known(:,1)-mean(val_known(:,1));
-    %iQ=inv(Q);
-    iQ=inv(Q+0.00000001*eye(size(Q,1))); % NEEDED FOR SOM UNSTABLE MATRIX INVERSION
+    iQ=inv(Q); % 
     try
         logdetQ=logdet(Q);
     catch
@@ -102,9 +91,10 @@
     sigma2_est=0;
 elseif method==3,
     mgstat_verbose(sprintf('%s : gauss likelihood',mfilename),-1);
-    %L=(-.5*dm'*inv(Cm)*dm);    
-    %logL = logdet(Cm)
-    L = 0.5*logdet(Cm) + (-.5*dm'*inv(Cm)*dm);
 
+    f1 = -.5*log(2*pi^nknown);
+    f2 = -0.5*logdet(Cm);
+    f3 = (-.5*dm'*inv(Cm)*dm); 
+    L = f1 + f2 +f3;
     sigma2_est=0;
 end
diff --git a/sgems/sgems_grid_py.m b/sgems/sgems_grid_py.m
@@ -204,16 +204,18 @@
     fname='ti.sgems';
     sgems_write_grid(1:1:ny,1:1:nx,1:1:nz,S.ti(:),fname,'ti','property');
     S.ti_file=fname;
-    S=rmfield(S,'ti');
+    S=rmfield(S,'ti');    
 end
 
+if isfield(S,'ti_file');
+    if exist(S.ti_file,'file')==2;
+        O_sgems=sgems_read(S.ti_file);
+        if ~isfield(S,'ti_property_id');S.ti_property_id=1;end
+        S.XML.parameters.PropertySelector_Training.grid=O_sgems.grid_name;
+        S.XML.parameters.PropertySelector_Training.property=O_sgems.property{S.ti_property_id};
+    end
+end
 
-if exist(S.ti_file,'file')==2;
-    O_sgems=sgems_read(S.ti_file);
-    if ~isfield(S,'ti_property_id');S.ti_property_id=1;end
-    S.XML.parameters.PropertySelector_Training.grid=O_sgems.grid_name;
-    S.XML.parameters.PropertySelector_Training.property=O_sgems.property{S.ti_property_id};
-end    
 
 
 % update marginal PDF for SNESIM if not set by user