compute.m

matlab中常用电法各种装置的电磁正演计算

function varargout = compute(varargin)

% compute.m is part of the CR1Dmod forward modeling package.
% 
% It contains functions relating to the compute GUI, from which calculation
% speciffic parameters are defined and the actual computation started. The 
% functions for storing and plotting results are also located in this file.
%
% Default values for the different parameters are located in the subfunction
% "initialize_gui", and may be modified by the user.
%
% Functions for plotting and storing the results are located at the end of
% the file.
%
% Written by:
% Thomas Ingeman-Nielsen
% The Arctic Technology Center, BYG
% Technical University of Denmark
% Email: [e-mail]

% Begin initialization code - DO NOT EDIT
gui_Singleton = 1;
gui_State = struct('gui_Name',       mfilename, ...
                   'gui_Singleton',  gui_Singleton, ...
                   'gui_OpeningFcn', @compute_OpeningFcn, ...
                   'gui_OutputFcn',  @compute_OutputFcn, ...
                   'gui_LayoutFcn',  [] , ...
                   'gui_Callback',   []);
if nargin & isstr(varargin{1})
    gui_State.gui_Callback = str2func(varargin{1});
end

if nargout
    [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});
else
    gui_mainfcn(gui_State, varargin{:});
end
% End initialization code - DO NOT EDIT

% --- Executes just before compute is made visible.
function compute_OpeningFcn(hObject, eventdata, handles, varargin)
% This function has no output args, see OutputFcn.
% hObject    handle to figure
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
% varargin   command line arguments to compute (see VARARGIN)


% Choose default command line output for compute
handles.output = hObject;

% Update handles structure
guidata(hObject, handles);

if ~isempty(varargin)
    if ishandle(varargin{1})
        handles.modelwin = varargin{1};  % handle to model window
    else
        return;          % bad input
    end
else
    return;              % bad input
end

guidata(hObject, handles);


if strcmp(get(hObject,'Visible'),'off')
    initialize_gui(hObject, handles);
end

handles = guidata(hObject);

setup_gui(hObject, [], handles, 'Default');
% UIWAIT makes compute wait for user response (see UIRESUME)
%uiwait(hObject);


% --- Outputs from this function are returned to the command line.
function varargout = compute_OutputFcn(hObject, eventdata, handles)
% varargout  cell array for returning output args (see VARARGOUT);
% hObject    handle to figure
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Get default command line output from handles structure

if isstruct(handles) & isfield(handles,'output')
    varargout{1} = handles.output;
else
    varargout{1} = [];
end

% -------------------------------------------------------------------------
function calculate_batch(hObject, eventdata, handles)

storeAllSettings(hObject, [], handles);

modelwin_han = guidata(handles.modelwin);

model.cparams = modelwin_han.cparams;
model.config = modelwin_han.config;
model.layers = modelwin_han.layers;

batmodels = modelwin_han.batch;

for k = 1:length(batmodels)
    Calculate_button_Callback(hObject, eventdata, handles, batmodels(k));
end
    
% -------------------------------------------------------------------------
function Calculate_button_Callback(hObject, eventdata, handles, varargin)

if nargin == 4 && isstruct(varargin{1})
    model = varargin{1};
    batch = 1;
else
    batch = 0;
    storeAllSettings(hObject, [], handles);
    
    modelwin_han = guidata(handles.modelwin);
    
    model.cparams = modelwin_han.cparams;
    model.config = rmfield(modelwin_han.config,{'plot_handle'});
    model.layers = rmfield(modelwin_han.layers,{'label_handle'});
    if isfield(model.layers,'interface_handle')
        model.layers = rmfield(model.layers,{'interface_handle',            ...
                'prop_lab_handle'});
    end
end

if ~isfield(model,'name')
    model.name = [];
end

for k = 1:length(model.layers)
    if isstr(model.layers(k).mu)
        model.layers(k).mu = str2double(model.layers(k).mu);
    end
end

if model.cparams.domain ~= 'DC'
    if isfield(model.cparams,'times')
        model.cparams.times = str2num(model.cparams.times);
    else
        model.cparams.freq = str2num(model.cparams.freq);
    end
end
result = [];

% prepare plotting
if ~isfield(handles,'fig') || isempty(handles.fig)       || ...
        ~get(handles.Reusefig_checkbox,'value')          || ...
        ~strcmp(model.config.type,handles.figtype{1})    || ...
        ~strcmp(model.cparams.domain,handles.figtype{2}) || ...
        ~ishandle(handles.fig)
    handles.fig = figure('visible','off');
    handles.figtype = {model.config.type;model.cparams.domain};
    guidata(hObject,handles);
else
    children = get(handles.fig,'children');
    ax = findobj(children, 'flat', 'type', 'axes');
    set(ax, 'nextplot', 'add');
end

switch model.config.type
    case {'TEM Central Loop'}
        result = temfwd(model.config, model.layers, model.cparams);
        if ~isempty(result)
            %plot_tem(result, model.config, model.cparams, handles.fig);
            plotdat('tem',result, model.config, model.cparams, handles.fig);
        end
    case {'HCP FDEM (HLEM)'}
        result = fdemfwd(model.config, model.layers, model.cparams);
        if ~isempty(result)
            %plot_hcp(result, model.config, model.cparams, handles.fig);
            plotdat('hcp',result, model.config, model.cparams, handles.fig);
        end % if
    otherwise
        switch model.cparams.domain
            case {'FD'}
                switch model.config.type
                    case {'Dipole-Dipole',                             ...
                                'General Surface Array'}
                        result = emgsafwd(model.config,                ...
                            model.layers, model.cparams, model.name);
                    case {'*Capacitance*', '*GSA capacitance*'}
                        result = emgsa_cap_fwd(model.config,                ...
                            model.layers, model.cparams);
                    otherwise
                        disp(['The configuration is not implemented'    ...
                                ' for frequency domain calculations.']);
                        return
                end
                if ~isempty(result)
                    if isfield(model, 'name')
                        [result.name] = deal(model.name);
                    end
                    %plot_emgsa(result, model.config,                    ...
                    %    model.cparams, handles.fig);    
                    plotdat('emgsa', result, model.config,               ...
                        model.cparams, handles.fig);    
                end % if
            case {'DC'}
                result = dcgsafwd(model.config, model.layers,          ...
                    model.cparams);
                if ~isempty(result)
                    plotdat('dcgsa', result, model.config, model.cparams, handles.fig);
                end % if
            otherwise
                disp(['Time domain calculations not yet implemented '   ...
                        'for this configuration']);        
        end % switch
end % switch

if ~isempty(result)
    if getappdata(0, 'debug')
        disp('compute.m:Calculate_button_Callback:debug:');
        disp('    Result variable exists and is not empty!');
    end
    
    if isfield(model.config,'plot_handle')
        rmfield(model.config,'plot_handle');
    end
    for k = 1:length(model.layers)
        if isfield(model.layers(k),'label_handle')
            rmfield(model.layers(k),'label_handle');
        end
        if isfield(model.layers(k),'interface_handle')
            rmfield(model.layers(k),'interface_handle');
        end
        if isfield(model.layers(k),'thick_lab_handle')
            rmfield(model.layers(k),'thick_lab_handle');
        end
    end
    model.result = result;
    % always save to 'forward.mat'
    save('forward.mat','model','-mat');

    % update plots before asking to save
    drawnow;
    if batch
        if ~isfield(model, 'name') || isempty(model.name)
            model.name = 'Model';
        end
        filename = [model.name '_' strrep(datestr(now),':','-')];
        filename = strrep(filename,' ','_');
        saveforward(model, filename);
    else
        saveforward(model);
    end
    %    For standalone mode
    %    assignin('base', 'model',model);
end


% -------------------------------------------------------------------------
function initialize_gui(hObject, handles)

handles.DefaultStrings = {                          ;                   ...
        handles.Times_edit,     'logspace(-6,-3,30)';                   ...
        handles.NST_tol_edit,   '1e-8'              ;                   ...
        handles.FST_err_edit,   '1e-8'              ;                   ...
        handles.Freq_edit,      'logspace(-1,4,25)' ;                   ...
        handles.NHT_tol_edit,   '1e-8'              ;                   ...
        handles.FHT_err_edit,   '1e-8'              ;                   ...
        handles.FDsp_NDEC_edit, '10'                ;                   ...
        handles.FDsp_Bmin_edit, '0.001'             ;                   ...
        handles.FDsp_Bmax_edit, '1000'              ;                   ...
        handles.Quad_tol_edit,  '1e-5'              ;                   ...
        handles.Rsp_NDEC_edit,  '10'                ;                   ...
        handles.Seg_tol_edit,   '1e-12'             ;                   ...
        handles.Max_seg_edit,   '300'               };

for k = 1:length(handles.DefaultStrings);
    set(cell2mat(handles.DefaultStrings(k,1)), 'String',                ...
        handles.DefaultStrings{k,2});
end % k

%PropName   =                      {'Enable', 'Value'};
Settings(1).Config  = {'TEM Central Loop'};
Settings(1).defaultDomain = 3; % = Time Domain
Settings(1).Defaults = {...        %  DC         FD	        TD     group
        handles.DC_radiobutton,      'off', [], 'off', 0 , 'off', 0 , 1 ;...  1
        handles.Freqdom_radiobutton, 'off', [], 'on' , 1 , 'on' , 0 , 1 ;...  2
        handles.Timedom_radiobutton, 'off', [], 'on' , 0 , 'on' , 1 , 1 ;...  3
        handles.Full_sol_radiobutton,'off', [], 'on' , 0 , 'on' , 0 , 2 ;...  4
        handles.Quasi_radiobutton,   'off', [], 'on' , 1 , 'on' , 1 , 2 ;...  5
        handles.Times_edit,          'off', [], 'off', [], 'on' , [], 3 ;...  6
        handles.Times_txt,           'off', [], 'off', [], 'on' , [], 3 ;...  7
        handles.NST_radiobutton,     'off', [], 'off', 0 , 'on' , 0 , 4 ;...  8
        handles.NST_tol_edit,        'off', [], 'off', [], 'off', [], 4 ;...  9
        handles.NST_tol_txt,         'off', [], 'off', [], 'off', [], 4 ;... 10
        handles.FST_radiobutton      'off', [], 'off', 0 , 'on' , 1 , 4 ;... 11
        handles.FST_err_edit,        'off', [], 'off', [], 'on' , [], 4 ;... 12