diff --git a/openair1/SIMULATION/TOOLS/random_channel.c b/openair1/SIMULATION/TOOLS/random_channel.c index 5c78a3a794..268d0166b1 100644 --- a/openair1/SIMULATION/TOOLS/random_channel.c +++ b/openair1/SIMULATION/TOOLS/random_channel.c @@ -128,6 +128,7 @@ void fill_channel_desc(channel_desc_t *chan_desc, chan_desc->first_run = 1; chan_desc->ip = 0.0; chan_desc->max_Doppler = max_Doppler; + chan_desc->Doppler_phase_cur = calloc(nb_rx, sizeof(double)); chan_desc->ch = calloc(nb_tx*nb_rx, sizeof(struct complexd *)); chan_desc->chF = calloc(nb_tx*nb_rx, sizeof(struct complexd *)); chan_desc->a = calloc(nb_taps, sizeof(struct complexd *)); @@ -1653,6 +1654,37 @@ channel_desc_t *new_channel_desc_scm(uint8_t nb_tx, 0); break; + case SAT_LEO_TRANS: + case SAT_LEO_REGEN: + nb_taps = 1; + Td = 0; + channel_length = 1; + ricean_factor = 0.0; + aoa = 0.0; + maxDoppler = 0; + chan_desc->sat_height = 600e3; + chan_desc->enable_dynamic_delay = true; + chan_desc->enable_dynamic_Doppler = false; // TODO: requires UE to support continuous Doppler estimation, compensation and pre-compensation + fill_channel_desc(chan_desc,nb_tx, + nb_rx, + nb_taps, + channel_length, + default_amp_lin, + NULL, + NULL, + Td, + sampling_rate, + channel_bandwidth, + ricean_factor, + aoa, + forgetting_factor, + maxDoppler, + channel_offset, + path_loss_dB, + 0); + printf("%s: satellite orbit height %f km\n", map_int_to_str(channelmod_names, channel_model), chan_desc->sat_height / 1000); + break; + default: LOG_W(OCM,"channel model not yet supported\n"); free(chan_desc); @@ -1711,6 +1743,7 @@ void free_channel_desc_scm(channel_desc_t *ch) { free(ch->R_sqrt[i]); free(ch->R_sqrt); + free(ch->Doppler_phase_cur); free(ch->ch); free(ch->chF); free(ch->a); @@ -1743,9 +1776,9 @@ int random_channel(channel_desc_t *desc, uint8_t abstraction_flag) { struct complexd phase, alpha, beta; start_meas(&desc->random_channel); - // For AWGN channel, the received signal (Srx) is equal to transmitted signal (Stx) plus noise (N), i.e., Srx = Stx + N, + // For AWGN and SAT_LEO_* channels, the received signal (Srx) is equal to transmitted signal (Stx) plus noise (N), i.e., Srx = Stx + N, // therefore, the channel matrix is the identity matrix. - if (desc->modelid == AWGN) { + if (desc->modelid == AWGN || desc->modelid == SAT_LEO_TRANS || desc->modelid == SAT_LEO_REGEN) { for (aarx=0; aarxnb_rx; aarx++) { for (aatx = 0; aatx < desc->nb_tx; aatx++) { desc->ch[aarx+(aatx*desc->nb_rx)][0].r = aarx%desc->nb_tx == aatx ? 1.0 : 0.0; @@ -2071,6 +2104,8 @@ static void display_channelmodel(channel_desc_t *cd,int debug, telnet_printfunc_ cd->forgetting_factor); prnt("Initial phase: %lf nb_path: %i \n", cd->ip, cd->nb_paths); + if (cd->modelid == SAT_LEO_TRANS || cd->modelid == SAT_LEO_REGEN) + prnt("satellite orbit height: %f\n", cd->sat_height); for (int i=0; inb_taps ; i++) { prnt("taps: %i lin. ampli. : %lf delay: %lf \n",i,cd->amps[i], cd->delays[i]); diff --git a/openair1/SIMULATION/TOOLS/sim.h b/openair1/SIMULATION/TOOLS/sim.h index 9f1e74dd1f..8297d9c365 100644 --- a/openair1/SIMULATION/TOOLS/sim.h +++ b/openair1/SIMULATION/TOOLS/sim.h @@ -126,6 +126,18 @@ typedef struct { char *model_name; /// flags to properly trigger memory free unsigned int free_flags; + /// time stamp when the time varying channel emulation starts (when client connected) + uint64_t start_TS; + /// height of LEO satellite + float sat_height; + /// flag to enable dynamic delay simulation for LEO satellite + bool enable_dynamic_delay; + /// flag to enable dynamic Doppler simulation for LEO satellite + bool enable_dynamic_Doppler; + /// Doppler phase increment (might vary over time, e.g. for LEO satellite) + float Doppler_phase_inc; + /// current Doppler phase of each RX antenna (for continuous phase from one block to the next) + float *Doppler_phase_cur; } channel_desc_t; typedef struct { @@ -218,6 +230,8 @@ typedef enum { EPA_low, EPA_medium, EPA_high, + SAT_LEO_TRANS, + SAT_LEO_REGEN, } SCM_t; #define CHANNELMOD_MAP_INIT \ {"custom",custom},\ @@ -253,6 +267,8 @@ typedef enum { {"EPA_low",EPA_low},\ {"EPA_medium",EPA_medium},\ {"EPA_high",EPA_high},\ + {"SAT_LEO_TRANS",SAT_LEO_TRANS},\ + {"SAT_LEO_REGEN",SAT_LEO_REGEN},\ {NULL, -1} #define CONFIG_HLP_SNR "Set average SNR in dB (for --siml1 option)\n" diff --git a/radio/rfsimulator/apply_channelmod.c b/radio/rfsimulator/apply_channelmod.c index 917a22924c..5bb207e5be 100644 --- a/radio/rfsimulator/apply_channelmod.c +++ b/radio/rfsimulator/apply_channelmod.c @@ -28,7 +28,7 @@ #include #include #include - +#include #include #include @@ -52,14 +52,85 @@ either we regenerate the channel (call again random_channel(desc,0)), or we keep it over subframes legacy: we regenerate each sub frame in UL, and each frame only in DL */ -void rxAddInput( const c16_t *input_sig, - c16_t *after_channel_sig, - int rxAnt, - channel_desc_t *channelDesc, - int nbSamples, - uint64_t TS, - uint32_t CirSize - ) { +void rxAddInput(const c16_t *input_sig, + c16_t *after_channel_sig, + int rxAnt, + channel_desc_t *channelDesc, + int nbSamples, + uint64_t TS, + uint32_t CirSize) +{ + if ((channelDesc->sat_height > 0) && (channelDesc->enable_dynamic_delay || channelDesc->enable_dynamic_Doppler)) { // model for transparent satellite on circular orbit + /* assumptions: + - The Earth is spherical, the ground station is static, and that the Earth does not rotate. + - An access or link is possible from the satellite to the ground station at all times. + - The ground station is located at the North Pole (positive Zaxis), and the satellite starts from the initial elevation angle 0° in the second quadrant of the YZplane. + - Satellite moves in the clockwise direction in its circular orbit. + */ + const double radius_earth = 6371e3; // m + const double radius_sat = radius_earth + channelDesc->sat_height; + const double GM_earth = 3.986e14; // m^3/s^2 + const double w_sat = sqrt(GM_earth / (radius_sat * radius_sat * radius_sat)); // rad/s + + // start_time and end_time are when the pos_sat_z == pos_gnb_z (elevation angle == 0 and 180 degree) + const double start_phase = -acos(radius_earth / radius_sat); // SAT is just rising above the horizon + const double end_phase = -start_phase; // SAT is just falling behind the horizon + const double start_time = start_phase / w_sat; // in seconds + const double end_time = end_phase / w_sat; // in seconds + const uint64_t duration_samples = (end_time - start_time) * channelDesc->sampling_rate; + const double t = start_time + ((TS - channelDesc->start_TS) % duration_samples) / channelDesc->sampling_rate; + + const double pos_sat_x = 0; + const double pos_sat_y = radius_sat * sin(w_sat * t); + const double pos_sat_z = radius_sat * cos(w_sat * t); + + const double vel_sat_x = 0; + const double vel_sat_y = w_sat * radius_sat * cos(w_sat * t); + const double vel_sat_z = -w_sat * radius_sat * sin(w_sat * t); + + const double pos_ue_x = 0; + const double pos_ue_y = 0; + const double pos_ue_z = radius_earth; + + const double dir_sat_ue_x = pos_ue_x - pos_sat_x; + const double dir_sat_ue_y = pos_ue_y - pos_sat_y; + const double dir_sat_ue_z = pos_ue_z - pos_sat_z; + + const double dist_sat_ue = sqrt(dir_sat_ue_x * dir_sat_ue_x + dir_sat_ue_y * dir_sat_ue_y + dir_sat_ue_z * dir_sat_ue_z); + const double vel_sat_ue = (vel_sat_x * dir_sat_ue_x + vel_sat_y * dir_sat_ue_y + vel_sat_z * dir_sat_ue_z) / dist_sat_ue; + + double dist_gnb_sat = 0; + double vel_gnb_sat = 0; + if (channelDesc->modelid == SAT_LEO_TRANS) { + const double pos_gnb_x = 0; + const double pos_gnb_y = 0; + const double pos_gnb_z = radius_earth; + + const double dir_gnb_sat_x = pos_sat_x - pos_gnb_x; + const double dir_gnb_sat_y = pos_sat_y - pos_gnb_y; + const double dir_gnb_sat_z = pos_sat_z - pos_gnb_z; + + dist_gnb_sat = sqrt(dir_gnb_sat_x * dir_gnb_sat_x + dir_gnb_sat_y * dir_gnb_sat_y + dir_gnb_sat_z * dir_gnb_sat_z); + vel_gnb_sat = (vel_sat_x * dir_gnb_sat_x + vel_sat_y * dir_gnb_sat_y + vel_sat_z * dir_gnb_sat_z) / dist_gnb_sat; + } + + const double c = 299792458; // m/s + const double prop_delay = (dist_gnb_sat + dist_sat_ue) / c; + if (channelDesc->enable_dynamic_delay) + channelDesc->channel_offset = prop_delay * channelDesc->sampling_rate; + + const double f_Doppler_shift_sat_ue = (vel_sat_ue / (c - vel_sat_ue)) * channelDesc->center_freq; + const double f_Doppler_shift_gnb_sat = (-vel_gnb_sat / c) * channelDesc->center_freq; + if (channelDesc->enable_dynamic_Doppler) + channelDesc->Doppler_phase_inc = 2 * M_PI * (f_Doppler_shift_gnb_sat + f_Doppler_shift_sat_ue) / channelDesc->sampling_rate; + + static uint64_t last_TS = 0; + if(TS - last_TS >= channelDesc->sampling_rate) { + last_TS = TS; + LOG_I(HW, "Satellite orbit: time %f s, Doppler: gNB->SAT %f kHz, SAT->UE %f kHz, Delay %f ms\n", t, f_Doppler_shift_gnb_sat / 1000, f_Doppler_shift_sat_ue / 1000, prop_delay * 1000); + } + } + // channelDesc->path_loss_dB should contain the total path gain // so, in actual RF: tx gain + path loss + rx gain (+antenna gain, ...) // UE and NB gain control to be added @@ -94,9 +165,16 @@ void rxAddInput( const c16_t *input_sig, } //l } - // Fixme: lround(), rount(), ... is detected by valgrind as error, not found why - out_ptr->r += lround(rx_tmp.r*pathLossLinear + noise_per_sample*gaussZiggurat(0.0,1.0)); - out_ptr->i += lround(rx_tmp.i*pathLossLinear + noise_per_sample*gaussZiggurat(0.0,1.0)); + if (channelDesc->Doppler_phase_inc != 0.0) { + double complex in = CMPLX(rx_tmp.r, rx_tmp.i); + double complex out = in * cexp(channelDesc->Doppler_phase_cur[rxAnt] * I); + rx_tmp.r = creal(out); + rx_tmp.i = cimag(out); + channelDesc->Doppler_phase_cur[rxAnt] += channelDesc->Doppler_phase_inc; + } + + out_ptr->r = lround(rx_tmp.r*pathLossLinear + noise_per_sample*gaussZiggurat(0.0,1.0)); + out_ptr->i = lround(rx_tmp.i*pathLossLinear + noise_per_sample*gaussZiggurat(0.0,1.0)); out_ptr++; } diff --git a/radio/rfsimulator/simulator.c b/radio/rfsimulator/simulator.c index 82561c9941..d944e028bb 100644 --- a/radio/rfsimulator/simulator.c +++ b/radio/rfsimulator/simulator.c @@ -709,6 +709,10 @@ static bool flushInput(rfsimulator_state_t *t, int timeout, int nsamps_for_initi samplesVoid[i]=(void *)&v; rfsimulator_write_internal(t, t->lastWroteTS > 1 ? t->lastWroteTS - 1 : 0, samplesVoid, 1, t->tx_num_channels, 1, false); + + buffer_t *b = &t->buf[conn_sock]; + if (b->channel_model) + b->channel_model->start_TS = t->lastWroteTS; } else { if ( events[nbEv].events & (EPOLLHUP | EPOLLERR | EPOLLRDHUP) ) { socketError(t,fd); @@ -761,6 +765,8 @@ static bool flushInput(rfsimulator_state_t *t, int timeout, int nsamps_for_initi nsamps_for_initial; LOG_D(HW, "UE got first timestamp: starting at %lu\n", t->nextRxTstamp); b->trashingPacket=true; + if (b->channel_model) + b->channel_model->start_TS = t->nextRxTstamp; } else if (b->lastReceivedTS < b->th.timestamp) { int nbAnt= b->th.nbAnt;