#!/usr/local/bin/rpl -sdp /* ================================================================================ Algorithme de l'Obèle Copyright 2001, BERTRAND Joël. ================================================================================ Entrées : néant -------------------------------------------------------------------------------- Sorties : néant -------------------------------------------------------------------------------- Effets de bord : néant ================================================================================ */ OBELE << rad 31 sf erase "" disp "Algorithme de l'obèle" disp "{\\Large\\sl Algorithme de l'obèle}" pr1 drop { "standard*(*)" } if "lambda" "existence" inquire then { { "name" "lambda" } "sequential" "replace" "writeonly" "formatted" } open else { { "name" "lambda" } "sequential" "new" "writeonly" "formatted" } open end format /* Paramètres d'entrée */ 4 // Nombre d'antennes 64 // Nombre de mobiles 64 // Facteur d'étalement true // Présence de bruit 1 // Seuil (contrainte à tenir en terme de C/I) 1E-5 // Niveau de bruit 1E-8 // Critère de convergence "Statistique" // Modèle de canal ("Statistique" ou "Aléatoire") .5 // Distance entre les différents capteurs de l'antenne { } // Directions des trajets (simple déclaration de variable) true // Initialisation omnidirectionnelle (ou par un contrôle // de puissance élémentaire) true // Normalisation des diagrammes true // Diagrammes d'antenne en coordonnées polaires true // Impression du résultat 3 // Nombre de paquets de mobiles (modèle Statistique) .02 // Dispersion en fraction de '2*PI' (modèle Statistique) true // Egalité des puissance émises (modèle Statistique) true // Un mobile par paquet tracé dans les diagrammes -> UNITE N_ANTENNES N_MOBILES FACTEUR_ETALEMENT ALGORITHME_BRUITE SEUIL BRUIT EPS MODELE_CANAL DIST DIRECTIONS INITIALISATION_OMNIDIRECTIONNELLE DIAGRAMME_NORMALISE DIAGRAMMES_POLAIRES AUTORISATION_IMPRESSION PAQUETS DISPERSION EQUIPUISSANCE TRACE_UN_MOBILE << "" disp "\\vskip 3ex\\noindent" pr1 drop "Configuration" pr1 drop "\\hrule\\vskip 1ex" pr1 drop cr "Nombre d'antennes : " N_ANTENNES ->str + pr1 disp cr "Nombre de mobiles : " N_MOBILES ->str + pr1 disp cr "Type de canal : " MODELE_CANAL ->str + pr1 disp PAQUETS 1 ->list 0 con -> REPARTITION << PAQUETS DISPERSION DIST EQUIPUISSANCE N_ANTENNES N_MOBILES MODELE_CANAL INITIALISATION_R 'DIRECTIONS' sto 'REPARTITION' sto if ALGORITHME_BRUITE then SEUIL CONVERSION_ALGORITHME_BRUITE end if INITIALISATION_OMNIDIRECTIONNELLE then N_ANTENNES 1 2 ->list 0 con { 1 1 } 1 put 1 N_MOBILES 1 - start dup next N_MOBILES ->list else dup FACTEUR_ETALEMENT ALGORITHME_BRUITE SEUIL OPTIMISATION_SIMPLE end N_MOBILES N_MOBILES 2 ->list 0 con -1 0 0 -> LISTE_R V_PONDERATION F AVP VP PUISSANCES_INITIALES << rclf 2 sci 1 N_MOBILES for J rclf std " Utilisateur " J ->str + disp stof LISTE_R J get disp "" disp next stof /* Calcul des rapports C/I initiaux */ LISTE_R V_PONDERATION 'F' { 1 1 } 1 N_MOBILES for J 1 N_MOBILES for K if J K same then 0 else V_PONDERATION K get dup trn LISTE_R J get rot * * { 1 1 } get end puti next next drop2 if ALGORITHME_BRUITE not then F regv max swap drop list-> drop -> M C << 1 N_MOBILES for J M J C 2 ->list get re next N_MOBILES 1 ->list ->array dup abs / >> else N_MOBILES idn N_MOBILES 1 2 ->list SEUIL BRUIT * con swap F - inv swap * re array-> list-> drop2 1 ->list ->array end 2 sci "\\vskip 3ex\\noindent" pr1 drop " Rapports C/I initiaux" pr1 disp "\\hrule\\vskip 1ex" pr1 drop "" disp dup array-> 1 get ->list 1 true -> AUTORISATION_CALCUL << do geti if 0 < then false 'AUTORISATION_CALCUL' sto end until dup 1 same end drop2 if AUTORISATION_CALCUL then V_PONDERATION over 720 DIAGRAMMES_POLAIRES TRACE_UN_MOBILE N_MOBILES N_ANTENNES DIST DIAGRAMME_NORMALISE DIRECTIONS PAQUETS REPARTITION DIAGRAMME if DIAGRAMMES_POLAIRES not then { "Azimut" "Puissance" } label end "Diagrammes avant optimisation" title persist prlcd cllcd dup -> PUISSANCES << 0 1 PUISSANCES size 1 get for I PUISSANCES I 1 ->list get + next >> 'PUISSANCES_INITIALES' sto if ALGORITHME_BRUITE not then 0 else BRUIT end N_MOBILES FACTEUR_ETALEMENT CALCUL_RAPPORTS_SIGNAUX_INTERFERENCE else 3 dropn " Résolution impossible du système" end pr1 disp "" disp >> /* Boucle principale */ "\\vskip 3ex\\noindent" pr1 drop " Minimisation de la plus grande valeur propre Lambda" disp " Minimisation de la plus grande valeur propre $\\lambda$" pr1 drop "\\hrule\\vskip 1ex" pr1 drop "" disp std while VP AVP - abs EPS > repeat /* Normalisation des pondérations */ 1 N_MOBILES for J LISTE_R J get V_PONDERATION J get FACTEUR_ETALEMENT ALGORITHME_BRUITE SEUIL NORMALISATION array-> drop N_ANTENNES 1 2 ->list ->array 'V_PONDERATION' swap J swap put next /* Calcul de la matrice F */ 'F' { 1 1 } 1 N_MOBILES for J 1 N_MOBILES for K if J K same then 0 else V_PONDERATION K get dup trn LISTE_R J get rot * * { 1 1 } get end puti next next drop2 /* Calcul du plus grand vecteur propre gauche de F */ F legv -> MATRICE << /* Réécriture de la fonction max pour éviter les */ /* erreurs numériques d'arrondis apparaissant avec */ /* les racines doubles du polynôme caractéristique */ do MATRICE max until over re 0 >= if dup not then swap 'MATRICE' swap 0 put swap drop end end >> list-> drop rot swap -> COLONNE << 1 N_MOBILES for J dup J COLONNE 2 ->list get re swap next drop >> VP 'AVP' sto N_MOBILES 1 2 ->list ->array swap re dup 'VP' sto pr1 UNITE over 1 ->list swap write "Lambda = " swap ->str + disp /* Normalisation du plus grand vecteur gauche de F */ dup abs / /* Calcul des matrices T */ -> PG << 1 N_MOBILES for J LISTE_R 1 get 0 con 1 N_MOBILES for K if J K same not then LISTE_R K get PG K 1 2 ->list get * + end next next >> N_MOBILES ->list /* Calcul des vecteurs propres généralisés des matrices T */ -> LISTE_T << 1 N_MOBILES for J LISTE_T J get LISTE_R J get gregv re min swap drop list-> drop -> COLONNE << 1 N_ANTENNES for K dup K COLONNE 2 ->list get swap next drop N_ANTENNES 1 ->list ->array dup abs / N_ANTENNES 1 2 ->list rdm 'V_PONDERATION' J rot put >> next >> end /* Normalisation des pondérations */ 1 N_MOBILES for J LISTE_R J get V_PONDERATION J get FACTEUR_ETALEMENT ALGORITHME_BRUITE SEUIL NORMALISATION array-> drop N_ANTENNES 1 2 ->list ->array 'V_PONDERATION' swap J swap put next /* Pondérations */ "" disp " Pondérations optimales" pr1 disp "\\hrule\\vskip 1ex" pr1 drop "" disp 1 N_MOBILES for J "W(" std J ->str + ") = " + 2 sci 'V_PONDERATION' J get N_ANTENNES 1 ->list rdm pr1 ->str + disp next /* Calcul des puissances par mobile nécessaires */ "" disp "\\vskip 3ex\\noindent" pr1 drop " Calcul des puissances par mobile nécessaires" pr1 disp "\\hrule\\vskip 1ex" pr1 drop "" disp if ALGORITHME_BRUITE not then F regv max swap drop list-> drop -> M C << 1 N_MOBILES for J M J C 2 ->list get re next N_MOBILES 1 ->list ->array dup abs / >> V_PONDERATION over 720 DIAGRAMMES_POLAIRES TRACE_UN_MOBILE N_MOBILES N_ANTENNES DIST DIAGRAMME_NORMALISE DIRECTIONS PAQUETS REPARTITION DIAGRAMME if DIAGRAMMES_POLAIRES not then { "Azimut" "Puissance" } label end "Diagrammes après optimisation" title persist prlcd else if VP 1 > then " Absence de solution physique !" else N_MOBILES idn N_MOBILES 1 2 ->list SEUIL BRUIT * con swap F - inv swap * re array-> list-> drop2 1 ->list ->array V_PONDERATION over 720 DIAGRAMMES_POLAIRES TRACE_UN_MOBILE N_MOBILES N_ANTENNES DIST DIAGRAMME_NORMALISE DIRECTIONS PAQUETS REPARTITION DIAGRAMME if DIAGRAMMES_POLAIRES not then { "Azimut" "Puissance" } label end "Diagrammes après optimisation" title persist prlcd end end V_PONDERATION over 2 ->list "resultat_obele" store { "graphique.eps" "postscript eps enhanced monochrome dashed" } lcd-> dup if dup type 2 same then pr1 else dup array-> 1 get ->list pr1 drop end if dup type 2 same not then "P = " swap ->str + disp "" disp LISTE_R V_PONDERATION rot dup -> PUISSANCES << 0 1 PUISSANCES size 1 get for I PUISSANCES I 1 ->list get + next >> if PUISSANCES_INITIALES 0 same not then rclf swap 3 fix PUISSANCES_INITIALES swap %ch neg ->str " %" + swap stof else drop "absurde" end "\\vskip 3ex\\noindent" pr1 drop " Rapports C/I finaux " "(amélioration de la puissance émise : " + swap ->str + ")" + pr1 disp if ALGORITHME_BRUITE not then 0 else BRUIT end N_MOBILES FACTEUR_ETALEMENT CALCUL_RAPPORTS_SIGNAUX_INTERFERENCE "\\hrule\\vskip 1ex" pr1 drop "" disp pr1 disp "" disp cllcd if AUTORISATION_IMPRESSION then print else erase end else drop disp "" disp erase cllcd end >> >> UNITE close >> " Temps CPU utilisé : " disp 2 fix time disp std >> /* ================================================================================ Calcul des rapports C/I pour chaque mobile ================================================================================ Entrées : 4: liste contenant les matrices R de chaque mobile 3: liste contenant les pondérations affectées à chaque mobile 2: vecteur contenant les puissances 1: sigma ** 2 -------------------------------------------------------------------------------- Sorties : 1: liste contenant les rapports C/I -------------------------------------------------------------------------------- Effets de bord : néant ================================================================================ */ CALCUL_RAPPORTS_SIGNAUX_INTERFERENCE << -> R W P SIGMA N_MOBILES FACTEUR_ETALEMENT << 1 N_MOBILES for I P I 1 ->list get W I get dup trn swap R I get swap * * * { 1 1 } get re SIGMA 1 N_MOBILES for J if I J same then cycle end P J 1 ->list get W J get dup trn swap R I get swap * * * { 1 1 } get re + next / FACTEUR_ETALEMENT * next N_MOBILES ->list >> >> /* ================================================================================ Fonction de normalisation des vecteurs W de telle sorte que trn(W)*R*W = 1 ================================================================================ Entrées : 2: matrice R 1: vecteur W -------------------------------------------------------------------------------- Sorties : 1: vecteur W normalisé -------------------------------------------------------------------------------- Effets de bord : néant ================================================================================ */ NORMALISATION << -> R W FACTEUR_ETALEMENT ALGORITHME_BRUITE SEUIL << W dup trn R FACTEUR_ETALEMENT * if ALGORITHME_BRUITE then SEUIL / end W * * abs sqrt / >> >> /* ================================================================================ Fonction renvoyant une liste contenant les différentes matrices R ================================================================================ Entrées : 3: nombre d'antennes (entier) 2: nombre de mobiles (entier) 1: nombre de trajets (entier) -------------------------------------------------------------------------------- Sorties : 2: liste contenant autant de matrices R qu'il y a de mobiles 1: directions des mobiles -------------------------------------------------------------------------------- Effets de bord : néant ================================================================================ */ INITIALISATION_R << "" disp " Initialisation des matrices d'autocorrélation du canal" disp "" disp "\\vskip 3ex\\noindent" pr1 drop "Positions et puissances des différents récepteurs" pr1 drop "\\hrule\\vskip 1ex" pr1 drop { } dup -> PAQUETS DISPERSION DIST EQUIPUISSANCE NA NM MODELE DIRECTIONS REPARTITION_INTERNE << rclf if MODELE "Statistique" same then PAQUETS 1 ->list 0 con 'REPARTITION_INTERNE' sto /* Génération de matrices de covariance du canal grâce au modèle du Statistique */ deg 4 fix { 3 9 } 0 con DIST 180 25 0 0 PAQUETS 1 ->list 0 con -> COEFF // Coefficients du modèle : // - ligne 1 : direction du trajet en degrés; // - ligne 2 : puissance du trajet en dB; // - ligne 3 : retard du trajet en ns. D // Distance entre deux capteurs consécutifs comptée en // longueur d'onde SECTEUR // Demi angle d'ouverture d'un secteur (en degrés) DTHETA // Paramètre d'ouverture du modèle (en degrés) GISMIN // Gisement le plus faible vu du secteur GISMAX // Gisement le plus grand vu du secteur ANGLES_MOYENS << SECTEUR neg DTHETA 4 * - 'GISMIN' sto SECTEUR DTHETA 2 * + 'GISMAX' sto 'COEFF' { 1 1 } 0 puti DTHETA 2 / puti DTHETA 2 / neg puti DTHETA 2 / 1 - puti 1 DTHETA 2 / - puti 2 DTHETA * puti -2 DTHETA * puti 3 DTHETA * puti 4 DTHETA * puti -2 puti -7 puti -7 puti -4 puti -4 puti -9 puti -10 puti -15 puti -20 puti 0 puti 0 puti 0 puti 310 puti 310 puti 710 puti 1090 puti 1730 puti 2510 puti drop2 0 -> CUMUL << 1 COEFF size 2 get for I 'COEFF' 2 I 2 ->list 'COEFF' over get 10 / alog dup 'CUMUL' sto+ put next 1 COEFF size 2 get for I 'COEFF' 2 I 2 ->list 'COEFF' over get CUMUL / put next >> // Calcul de la répartition des mobiles dans les paquets 'REPARTITION_INTERNE' { 1 } 1 PAQUETS for P rand puti next drop2 REPARTITION_INTERNE array-> 1 get 2 swap for P + next -> CLEF << 1 PAQUETS for P 'REPARTITION_INTERNE' dup P 1 ->list get NM * CLEF / ip P 1 ->list swap if dup 1 < then drop 1 end put next >> 0 1 PAQUETS for P 'REPARTITION_INTERNE' P 1 ->list get + next NM - -> DIFFERENCE << if DIFFERENCE 0 > then while DIFFERENCE repeat rand PAQUETS * ip 1 + 1 ->list dup if 'REPARTITION_INTERNE' swap get dup 1 > then 1 - 'REPARTITION_INTERNE' rot rot put 'DIFFERENCE' 1 sto- else drop2 end end else while DIFFERENCE repeat rand PAQUETS * ip 1 + 1 ->list dup if 'REPARTITION_INTERNE' swap get dup NM < then 1 + 'REPARTITION_INTERNE' rot rot put 'DIFFERENCE' 1 sto+ else drop2 end end end >> 'ANGLES_MOYENS' { 1 } 1 PAQUETS for P rand 2 SECTEUR * * SECTEUR - puti next drop2 rclf std "Répartition : " REPARTITION_INTERNE ->str + pr1 "\\hrule\\vskip 1ex" pr1 drop disp "" disp stof // Boucle sur les paquets de mobiles 1 PAQUETS for P // Boucle sur les mobiles 1 REPARTITION_INTERNE P 1 ->list get for K rand 2 SECTEUR * * SECTEUR - DISPERSION * ANGLES_MOYENS P 1 ->list get + DIRECTIONS over 1 ->list + 'DIRECTIONS' sto "Azimut : " over ->hms ->str + dup " ° (HMS)" + disp "\\degre (HMS)" + cr pr1 drop NA COEFF size 2 get 2 ->list 0 con COEFF size 2 get dup 2 ->list 0 con 0 -> AZIMUT MD P TRAJETS_RETENUS << 1 COEFF size 2 get for I 'COEFF' 1 I 2 ->list get AZIMUT + -> G << if G SECTEUR <= G SECTEUR neg >= and SECTEUR 180 >= or then 1 NA for J 'MD' 2 i pi * * ->num D J 1 - * * G sin * exp J I 2 ->list swap put next 'P' COEFF 2 I 2 ->list get I dup 2 ->list swap put 1 'TRAJETS_RETENUS' sto+ end >> next MD P if EQUIPUISSANCE then 1 else nrand sq end "Puissance : " over ->str + cr pr1 disp * over trn * * // Rajout de bruit pour éviter d'avoir une // matrice R de rang non plein if TRAJETS_RETENUS over size 1 get < then dup idn over abs 1E-6 / * + end >> next next >> rad else /* Génération de matrices R aléatoires */ // Nombre de trajets 4 -> NT << // Boucle sur les mobiles 1 NM for K " Utilisateur " std K ->str + disp "" disp NT NA 2 ->list 0 con // Boucle sur les trajets 1 NT for L L 1 2 ->list rand 2 pi ->num * * DIRECTIONS over r->d 1 ->list + 'DIRECTIONS' sto nrand sq -> D P << std " Trajet " L ->str + disp 4 sci " -> Puissance " P ->str + disp 4 fix " -> Azimut " D r->d ->hms ->str + "° (HMS)" + disp "" disp // Boucle sur les capteurs 1 NA for J i ->num D sin J 1 - * * DIST * 2 pi ->num * 2E9 * * exp P * puti next drop >> next trn conj dup trn * next >> NM 1 ->list 1 con 'REPARTITION_INTERNE' sto end NM ->list swap stof REPARTITION_INTERNE DIRECTIONS >> "" disp >> /* ================================================================================ Fonction permettant de convertir l'algorithme non bruité en sa version bruitée ================================================================================ Entrées : 2: liste contenant les différentes matrices R 1: seuil -------------------------------------------------------------------------------- Sorties : 1: liste contenant les nouvelles matrices R' -------------------------------------------------------------------------------- Effets de bord : néant ================================================================================ */ CONVERSION_ALGORITHME_BRUITE << -> L_R SEUIL << 1 L_R size for J 'L_R' dup J get SEUIL * J swap put next L_R >> >> /* ================================================================================ Calcul simple des pondérations et des puissances ================================================================================ Entrées : 1: liste contenant les différentes matrices R -------------------------------------------------------------------------------- Sorties : 2: liste contenant les pondérations et les puissances 1: valeur propre -------------------------------------------------------------------------------- Effets de bord : néant ================================================================================ */ OPTIMISATION_SIMPLE << -> LISTE FACTEUR_ETALEMENT ALGORITHME_BRUITE SEUIL << 1 LISTE size for N LISTE N get regv max swap drop list-> drop over size 1 get 1 2 ->list 0 con -> INDICE TABLEAU << 1 over size 1 get for M dup M 1 2 ->list get 'TABLEAU' swap M 1 2 ->list swap put next drop LISTE N get TABLEAU FACTEUR_ETALEMENT ALGORITHME_BRUITE SEUIL NORMALISATION >> next LISTE size ->list >> >> /* ================================================================================ Calcul du diagramme de rayonnement du réseau d'antennes ================================================================================ Entrées : 3: liste contenant tous les vecteurs de pondération 2: puissances 1: nombre de points à calculer par diagramme -------------------------------------------------------------------------------- Sorties : néant -------------------------------------------------------------------------------- Effets de bord : néant ================================================================================ */ DIAGRAMME << 0 -> PONDERATIONS PUISSANCES NB_POINTS DIAGRAMMES_POLAIRES TRACE_UN_MOBILE N_MOBILES N_ANTENNES DIST DIAGRAMME_NORMALISE DIRECTIONS PAQUETS REPARTITION MAXIMUM << cllcd { { 60 "ticsonly" 2 } { "automatic" "ticsonly" 10 } } axes if DIAGRAMMES_POLAIRES then 1 d->r -> PAS << 0 1 N_MOBILES for I 0 2 pi ->num * for T if DIST N_ANTENNES T PONDERATIONS I DIAGRAMME_NORMALISE PUISSANCES FCT_DIAGRAMME dup 3 pick > then swap end drop PAS step next >> dup 'MAXIMUM' sto dup r->c dup pmax neg pmin parametric { T 0 'MAXIMUM' } indep MAXIMUM res << T I DIRECTIONS FCT_DIRECTIONS >> steq 1 N_MOBILES for I draw next polar { T 0 '2*PI' } indep 2 pi ->num * NB_POINTS / res << DIST N_ANTENNES T PONDERATIONS I DIAGRAMME_NORMALISE PUISSANCES FCT_DIAGRAMME >> steq if TRACE_UN_MOBILE then 1 -> I << 1 PAQUETS for POINTEUR draw 'I' REPARTITION POINTEUR 1 ->list get ip sto+ next >> else 1 N_MOBILES for I draw next end else { X Y } autoscale 'Y' logscale parametric { T '-PI' 'PI' } indep 2 pi ->num * NB_POINTS / res << T r->d DIST N_ANTENNES T PONDERATIONS I DIAGRAMME_NORMALISE PUISSANCES FCT_DIAGRAMME r->c >> steq if TRACE_UN_MOBILE then 1 -> I << 1 PAQUETS for POINTEUR draw 'I' REPARTITION POINTEUR 1 ->list get ip sto+ next >> else 1 N_MOBILES for I draw next end end >> drax >> FCT_DIRECTIONS << -> T I DIRECTIONS << rclf deg DIRECTIONS I get dup cos T * swap sin T * i ->num * + swap stof >> >> FCT_DIAGRAMME << -> DIST N_ANTENNES T PONDERATIONS I DIAGRAMME_NORMALISE PUISSANCES << N_ANTENNES 1 ->list i ->num 2 pi ->num * * T sin DIST * * con 1 N_ANTENNES for J dup J 1 ->list get J 1 - * exp J 1 ->list swap put next PONDERATIONS I get array-> 1 get 1 ->list ->array swap dot abs 2 / if DIAGRAMME_NORMALISE not then PUISSANCES I 1 ->list get * end >> >> // vim: ts=4