Home > atmlab > randomize > iaaft > iaaft_loop_2d_vertical.m

iaaft_loop_2d_vertical

PURPOSE ^

This function is the main loop of the Iterative Amplitude Adapted Fourier

SYNOPSIS ^

function [surrogate, error_amplitude, error_spec] = iaaft_loop_2d_vertical(coeff_2d, sorted_values_prof)

DESCRIPTION ^

 This function is the main loop of the Iterative Amplitude Adapted Fourier
 Transform  (IAAFT) method to make surrogate fields. Use the m-files that
 start with surrogate* to run this function.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

DOWNLOAD ^

iaaft_loop_2d_vertical.m

SOURCE CODE ^

0001 function [surrogate, error_amplitude, error_spec] = iaaft_loop_2d_vertical(coeff_2d, sorted_values_prof)
0002 
0003 % This function is the main loop of the Iterative Amplitude Adapted Fourier
0004 % Transform  (IAAFT) method to make surrogate fields. Use the m-files that
0005 % start with surrogate* to run this function.
0006 
0007 % This Matlab version was written by Victor Venema,
0008 % Victor.Venema@uni-bonn.de, http:\\www.meteo.uni-bonn.de\victor, for the
0009 % generation of surrogate cloud fields. First version: May 2003.
0010 
0011 % INPUT:
0012 % coeff_2d:           The 2 dimensional Fourier coefficients that describe the
0013 %                     structure and implicitely pass the size of the matrix
0014 % sorted_values_prof: A vector with the wanted amplitudes (e.g. LWC of LWP
0015 %                     values) sorted in acending order for each height level.
0016 % OUTPUT:
0017 % y:                  The 2D IAAFT surrogate field
0018 % error_amplitude:    The amount of addaption that was made in the last
0019 %                     amplitude addaption relative to the total standard deviation.
0020 % error_spec:         The amont of addaption that was made in the last fourier
0021 %                     coefficient addaption relative to the total standard deviation
0022 
0023 % settings
0024 errorThresshold = 0.0012; % The addaptation made in the last iterative step relative to the total standard deviation
0025 timeThresshold  = Inf;    % 180;    % CPU-time in seconds that the iteration maximally uses
0026 speedThresshold = 1e-4    % 1e-8;   % Minimal convergence speed in the maximum error.
0027 
0028 % Initialse function
0029 [no_values_x, no_values_z] = size(coeff_2d);
0030 error_amplitude = 1;
0031 error_spec      = 1;
0032 oldTotalError   = 10;
0033 speed = 1;
0034 standard_deviation = 0;
0035 standard_deviation = std(sorted_values_prof(:));
0036 t=cputime;
0037 
0038 % The method starts with a randomized uncorrelated time series y with the saame pdf of
0039 % sorted values at each height level.
0040 surrogate = zeros(no_values_x, no_values_z);
0041 for j = 1:no_values_z
0042     [dummy, index] = sort(rand(size(sorted_values_prof(:,j))));
0043     surrogate(index,j) = sorted_values_prof(:,j);
0044 end
0045 surrogate = reshape(surrogate, no_values_x, no_values_z);
0046 
0047 % Main intative loop
0048 while ( (error_amplitude > errorThresshold | error_spec > errorThresshold) & (cputime-t < timeThresshold) & (speed > speedThresshold) ) %0.0001 300
0049     % addapt the power spectrum
0050     old_surrogate = surrogate;    
0051     x = ifft2(double(surrogate));
0052     phase = angle(x);
0053     x = double(coeff_2d) .* exp(i*phase);
0054     surrogate = fft2(double(x));
0055     
0056     difference = mean(mean(mean(abs(double(real(surrogate)) - double(real(old_surrogate))))));
0057     error_spec = difference/standard_deviation
0058       
0059     % addept the amplitude distribution
0060     old_surrogate = surrogate;
0061     surrogate = reshape(surrogate, no_values_x, no_values_z);
0062     for j = 1:no_values_z
0063         [dummy, index]     = sort(real(surrogate(:,j)));
0064         surrogate(index,j) = sorted_values_prof(:,j);
0065     end
0066     surrogate  = reshape(surrogate, no_values_x, no_values_z);
0067     
0068     difference = mean(mean(mean(abs(double(real(surrogate)) - double(real(old_surrogate))))));
0069     error_amplitude = difference / standard_deviation
0070 
0071     totalError = error_spec + error_amplitude
0072     speed = (oldTotalError - totalError) / totalError
0073     oldTotalError = totalError;
0074 end
0075 
0076 error_spec
0077 error_amplitude
0078 surrogate = real(surrogate);

Generated on Mon 15-Sep-2014 13:31:28 by m2html © 2005