Code covered by the BSD License

# The Gravity Perturbed Hohmann Transfer

### David Eagle (view profile)

28 Feb 2013 (Updated )

MATLAB script for solving the Hohmann transfer problem perturbed by non-spherical Earth gravity.

[fid, alt1, alt2, inc1, inc2, raan1, jdutc0, gst0, ...
```function [fid, alt1, alt2, inc1, inc2, raan1, jdutc0, gst0, ...
req, mu, lgrav mgrav, gravfname] = readdata(filename)

% read orbital elements and simulation
% control data file

% required by phohmann.m

% input

%  filename = name of oota data file

% output

%  fid = file id

%  oev1 = initial orbit orbital elements array
%  oev2 = final orbit orbital elements array

% Orbital Mechanics with MATLAB

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% open data file

fid = fopen(filename, 'r');

% check for file open error

if (fid == -1)

clc; home;

fprintf('\n\n  error: cannot find this file!!');

keycheck;

return;

end

% read 67 lines of data

for i = 1:1:67

cline = fgetl(fid);

switch i

case 10

% initial UTC calendar date

tl = size(cline);

ci = findstr(cline, ',');

% extract month, day and year

month = str2double(cline(1:ci(1)-1));

day = str2double(cline(ci(1)+1:ci(2)-1));

year = str2double(cline(ci(2)+1:tl(2)));

case 15

% initial UTC epoch

utstr = cline;

tl = size(utstr);

ci = findstr(utstr, ',');

% extract hours, minutes and seconds

utc_hr = str2num(utstr(1:ci(1)-1));

utc_min = str2num(utstr(ci(1)+1:ci(2)-1));

utc_sec = str2num(utstr(ci(2)+1:tl(2)));

% compute utc julian date

day = day + utc_hr / 24.0 + utc_min / 1440.0 + utc_sec / 86400.0;

jdutc0 = julian(month, day, year);

gst0 = gast1 (jdutc0);

case 23

alt1 = str2double(cline);

case 27

inc1 = str2double(cline);

case 31

raan1 = str2double(cline);

case 39

alt2 = str2double(cline);

case 43

inc2 = str2double(cline);

case 51

mu = str2double(cline);

case 55

req = str2double(cline);

case 59

gravfname = cline;

case 63

lgrav = str2double(cline);

case 67

mgrav = str2double(cline);

end
end

status = fclose(fid);

```