function [T,X] = air_mass(trange,N)

t = linspace(trange(1), trange(2), N);
% alpha = 1.3e-6; beta = 5.95e-5;
% Do2 = 3.36; POi = 98;
M0 = [140, 142, 0.0123, 0.0123]; 
% R_R = 1.17; Q_v = 2.97;
% RFt = 0.0154; V0 = 350;
% 
% Q_air = Q_v*sin(2*pi*R_R*t);
% qi = RFt*Q_air;
% A = qi/(2*pi*R_R);
% V_alv = A.*(1-cos(2*pi*R_R*t))+ V0;
% 
% P_alv = (beta*M)./V_alv;

% m = (q/beta)*P_alv - alpha*Do2*(P_alv-POi);

options = odeset('RelTol', 1e-08, 'AbsTol', 1e-08, 'MaxStep', 0.01);
tspan = t;

[T,X]=ode15s(@mass, tspan, M0, options);
% [~,~,P_alv,V_alv,Q_air] = mass(T,m);