computemripeaks.m

低矮房屋风压系数、风荷载计算分析matlab程序

function MRI_struct = computeMRIpeaks( bldg_struct, DIF_struct_cell, HURR_struct, MRI_list )

% check to make sure that bldg_struct and DIF_struct correspond
% if the terrain conditions are directional, then DIF_struct should be
% a structure array of length 2, with one "open" and one "suburban" 
% DIF_struct should be consistent


% PROCEDURE:


% for each building orientation:

% define list of wind directions relative to reference axis

% compute DIFs for each wind direction

err_dlgname = 'Error computing MRI peaks';

tol = 10^-3;
% Convert the input 'DIF_struct_cell' to a structure array 'DIF_struct' in which each element has the same fields:
if length(DIF_struct_cell)==1
    DIF_struct = DIF_struct_cell{1};
    n_r = length(DIF_struct(1).resp_names); % number of responses
    % make sure terrain in DIF results matches building input file:
    if ~ismember(bldg_struct.terrain, {DIF_struct.terrain,'Directional'})
        err = errordlg(['The terrain in the building input file (''' bldg_struct.terrain ...
            ''') does not match the terrain in the DIF results file (''' DIF_struct.terrain ''').'],err_dlgname);
        uiwait(err); MRI_struct.error = 1; return;
    end
elseif length(DIF_struct_cell)==2
    % Make sure building input file specifies 'Directional' terrain:
    if ~strcmp(bldg_struct.terrain,'Directional')
        err = errordlg(['Only one DIF results file expected for the terrain specified in the building input file (''' bldg_struct.terrain ...
            '''). Multiple DIF results allowed for ''Directional'' terrain only.' ],err_dlgname);
        uiwait(err); MRI_struct.error = 1; return;
    end
    fields1 = fieldnames(DIF_struct_cell{1});
    fields2 = fieldnames(DIF_struct_cell{2});
    fields = intersect(fields1,fields2);
    % Create a structure array containing only the common fields:
    for i=1:2
        for j=1:length(fields)
            DIF_struct(i).(fields{j}) = DIF_struct_cell{i}.(fields{j});
        end
    end
    % Make sure the terrain does not coincide for the two DIF results (one should be 'Open_Country' and one should be 'Suburban'):
    if strcmp(DIF_struct(1).terrain, DIF_struct(2).terrain)
        err = errordlg(['Both DIF files correspond to the same type of terrain (''' DIF_struct(1).terrain ...
            '''). DIF files for both ''Open_Country'' and ''Suburban'' terrain are required to compute MRIs by interpolation.'],err_dlgname);
        uiwait(err); MRI_struct.error = 1; return;
    % make sure units, wind directions, frame locations, and building dimensions coincide:
    elseif ~strcmp(DIF_struct(1).force_units,DIF_struct(2).force_units)
        err = errordlg('The force units must coincide for all DIF results used in interpolation', err_dlgname);
        uiwait(err); MRI_struct.error = 1; return;
    elseif ~strcmp(DIF_struct(1).length_units,DIF_struct(2).length_units)
        err = errordlg('The length units must coincide for all DIF results used in interpolation', err_dlgname);
        uiwait(err); MRI_struct.error = 1; return;
    elseif ~strcmp(DIF_struct(1).ws_units,DIF_struct(2).ws_units)
        err = errordlg('The wind speed units must coincide for all DIF results used in interpolation', err_dlgname);
        uiwait(err); MRI_struct.error = 1; return;
    elseif length(DIF_struct(1).theta)~=length(DIF_struct(i).theta) || any(DIF_struct(1).theta~=DIF_struct(i).theta)
        err = errordlg( 'The wind directions must coincide for all DIF results used in interpolation', err_dlgname);
        uiwait(err); MRI_struct.error = 1; return;
    elseif any(size(DIF_struct(1).frame_coords)~=size(DIF_struct(2).frame_coords)) || max(max(abs(DIF_struct(1).frame_coords-DIF_struct(2).frame_coords)))>tol
        err = errordlg( 'The frame coordinates must coincide for all DIF results used in interpolation', err_dlgname);
        uiwait(err); MRI_struct.error = 1; return;
    elseif max(max(abs(DIF_struct(1).d0-DIF_struct(2).d0)))>tol
        err = errordlg( 'The building dimensions must coincide for all DIF results used in interpolation',err_dlgname);
        uiwait(err); MRI_struct.error = 1; return;
    elseif size(DIF_struct(1).d,1)~=size(DIF_struct(2).d,1) || any(any(DIF_struct(1).d~=DIF_struct(2).d))
        err = errordlg( 'The model dimensions must coincide for all DIF results used in interpolation.',err_dlgname);
        uiwait(err); MRI_struct.error = 1; return;
    end
    % make sure the number of responses, the response names, and the response units coincide:
    if length(DIF_struct(1).resp_names)~=length(DIF_struct(2).resp_names)
        err = errordlg('The number of response quantities must coincide for all DIF results used in interpolation', err_dlgname);
        uiwait(err); MRI_struct.error = 1; return;
    end
    n_r = length(DIF_struct(1).resp_names);
    for j=1:n_r
        if ~strcmp(DIF_struct(1).resp_names(j), DIF_struct(2).resp_names(j))
            err = errordlg('The response names must coincide for all DIF results used in interpolation', err_dlgname);
            uiwait(err); MRI_struct.error = 1; return;
        elseif ~strcmp(DIF_struct(1).resp_units(j), DIF_struct(2).resp_units(j))
            err = errordlg('The response units must coincide for all DIF results used in interpolation', err_dlgname);
            uiwait(err); MRI_struct.error = 1; return;
        end
    end
else
    err = errordlg(['Unexpected number of DIF files (' num2str(length(DIF_struct_cell)) '; 1 or 2 expected.'], err_dlgname);
    uiwait(err); MRI_struct.error = 1; return;
end
MRI_struct.length_units = bldg_struct.length_units;
MRI_struct.force_units = bldg_struct.force_units;
MRI_struct.ws_units = bldg_struct.ws_units;
MRI_struct.d0 = bldg_struct.d0;
MRI_struct.d = DIF_struct(1).d;
MRI_struct.frame_coords = bldg_struct.frame_coords;
MRI_struct.terrain = bldg_struct.terrain;
MRI_struct.z0 = bldg_struct.z0;
MRI_struct.model_terrain = {DIF_struct.terrain};
MRI_struct.resp_names = DIF_struct(1).resp_names;
MRI_struct.resp_units = DIF_struct(1).resp_units;

% Identify which quantities are available to interpolate:
% initialize boolean values for which quantities to interpolate:
max_obs = 1;
min_obs = 1;
max_est = 1;
min_est = 1;

for i = 1:length(DIF_struct)
    if ~isfield(DIF_struct, 'X_max_obs') || isempty(DIF_struct(i).X_max_obs)
        max_obs = 0;
    elseif max_obs==1
        X_max_obs = DIF_struct(i).X_max_obs;
        n_f = size(X_max_obs,3);
        % 'wrap' DIFs above 360 and below 0 to facilitate resampling:
        DIF(i).X_max_obs_wrap = cat(2,X_max_obs,X_max_obs,X_max_obs);
    end

    if ~isfield(DIF_struct, 'X_min_obs') || isempty(DIF_struct(i).X_min_obs)
        min_obs = 0;
    elseif min_obs==1
        X_min_obs = DIF_struct(i).X_min_obs;
        n_f = size(X_min_obs,3);
        % 'wrap' DIFs above 360 and below 0 to facilitate resampling:
        DIF(i).X_min_obs_wrap = cat(2,X_min_obs,X_min_obs,X_min_obs);
    end

    if ~isfield(DIF_struct, 'X_max_est') || isempty(DIF_struct(i).X_max_est)
        max_est = 0;
    elseif max_est==1
        X_max_est = DIF_struct(i).X_max_est;
        n_f = size(X_max_est,3);
        % 'wrap' DIFs above 360 and below 0 to facilitate resampling:
        DIF(i).X_max_est_wrap = cat(2,X_max_est,X_max_est,X_max_est);
    end

    if ~isfield(DIF_struct, 'X_min_est') || isempty(DIF_struct(i).X_min_est)
        min_est = 0;
    elseif min_est==1
        X_min_est = DIF_struct(i).X_min_est;
        n_f = size(X_min_est,3);
        % 'wrap' DIFs above 360 and below 0 to facilitate resampling:
        DIF(i).X_min_est_wrap = cat(2,X_min_est,X_min_est,X_min_est);
    end
end
if ~max_obs & ~min_obs & ~max_est & ~min_est
    err = errordlg('Neither observed nor estimated DIF results are available for interpolation',err_dlgname);
    uiwait(err); MRI_struct.error = 1; return;
end

MRI_struct.max_obs = max_obs;
MRI_struct.min_obs = min_obs;
MRI_struct.max_est = max_est;
MRI_struct.min_est = min_est;

MRI_struct.DIF_struct = DIF_struct;

theta = DIF_struct(1).theta;

if any(theta~=sort(theta))
    error('DIF_struct.theta must be sorted in increasing order');
end    
% 'wrap' wind directions above 360 and below 0 to facilitate resampling:
theta_wrap = [theta(:)-360; theta(:); theta(:)+360];

% define upper and lower boundaries of sector for each wind direction:
theta_wrap_max = theta_wrap; theta_wrap_max(1:end-1) = theta_wrap(1:end-1)+ 0.5*diff(theta_wrap);
theta_wrap_min = theta_wrap; theta_wrap_min(2:end)   = theta_wrap(2:end)  - 0.5*diff(theta_wrap);

% define list of orientations of ridge axis of building:
d_a0 = 15; % increment in orientation (degrees)
alpha_0 = 0:d_a0:360-d_a0; % orientations (degrees clockwise from north)
n_o = length(alpha_0); % number of orientations considered

% Scale directional wind speeds to eave height over appropriate terrain: