Python3, 479 bytes

E=enumerate
def f(b):
 d,p={},[]
 for x,r in E(b):
  for y,v in E(r):d[v]=d.get(v,[])+[(x,y)]
 X,Y=d['O'][0]
 for x,y in d['#']:
  J,K=abs(j:=x-X),abs(k:=y-Y)
  if j and k:
   if J==K:
    if x<X:p+=[[(X,y),(x,Y)][y<Y]]
    else:p+=[[(x,Y),(X,y)][y<Y]]
   else:p+=[[(X,y),(x,Y)][J>K]]
  else:p+=[(x,y)]
 a,b,c,d=' ',' ',' ',' '
 for x,y in p:
  if y==Y and x<X:a='#'
  if y==Y and x>X:b='#'
  if x==X and y<Y:c='#'
  if x==X and y>Y:d='#'
 return f' {a} \n'+f'{c}O{d}\n'+f' {b} '

Try it online!

Perl, 156 153 148

$n=$_+1;
s/.+\n//;
$o=index$_,O;
$x=$o%$n-$-[0]%$n,
pos=$o+(($==($o-$-[0])/$n)&&$x/$=>-1&&$x/$=<=1?$=>0?-$n:$n:0<=>$x),
s/\G./x/ while s/#/./;
y/x/#/

Newlines for readability. Run with -0777p (5 added to count), input through STDIN.

More readable:

$n=$_+1;
s/.+\n//;
$o=index$_,O;
while(s/#/./){
    $x=$o%$n-$-[0]%$n;
    $==($o-$-[0])/$n;
    pos=$o+(
        $=
        &&$x/$=>-1
        &&$x/$=<=1
            ?$=>0
                ?-$n
                :$n
            :0<=>$x

    );
    s/\G./x/
}
y/x/#/

Mathematica, 206 bytes

f=(m=Characters@StringSplit@#;o=#-Join@@m~(p=Position)~"O"&/@m~p~"#";h=If[MemberQ[#3/@o,{m_,n_}/;m#<0&&m#<=n#2<-m#],"#",e]&;{e=".",h[1,1,#&],e,n="\n",h[1,-1,r=Reverse],"O",h[-1,1,r],n,e,h[-1,-1,#&],e}<>"")&

Somewhat more readable:

f[input_] := (
  map = Characters @ StringSplit @ input;
  planet = Flatten[Position[map, "O"]];
  objects = Map[# - planet &, Position[map, "#"]];
  helper[pattern_] := If[
    Length[Position[objects, pattern]] > 0, 
    "#", 
    "."
  ];
  StringJoin[{
    ".", h[{m_, n_} /; m < 0 && m <= n < -m], ".", "\n",
    h[{m_, n_} /; n < 0 && n < m <= -n], "O", h[{m_, n_} /; n > 0 && -n <= m < n], "\n",
    ".", h[{m_, n_} /; m > 0 && -m < n <= m], "."
  }]
);

As you can see, I only take a single input which is the string containing the input map. This is because I wouldn't use the size integer anyway. If you insist on the input format, I can change the input string to the second argument for one more character.

As for the algorithm, I'm just collecting all the {Δx, Δy} pairs of the objects and look for one pair in each quadrant with the conditions passed to h. I suspect this is also the part that could be golfed down further.

perl, 270 chars

3 steps: 1) find the planet and satellites in the input, 2) determine quadrant for each satellite, 3) output

sub _{"#"}$D=<>;$\=$/;for$L(1..$D){for$M(1..$D){$/=\1;$_=<>;$p=$M,$q=$L if/O/;push@s,[$M,$L] if/#/}$/=$\;<>}for(@s){$x=$_->[0]-$p,$y=$_->[1]-$q;($x<0&&$x<$y&&$y<=-$x?$l:$x>0&&-$x<=$y&&$y<$x?$r:$y<0&&$y<=$x&&$x<-$y?$u:$d)=_}$_=$l?" #":_;$u&&print;print $l."O$r";$d&&print

Less golfed:

sub _{"#"}
$D=<>;
$\=$/;
for$L(1..$D){
  $/=\1;
  for$M(1..$D){
    $_=<>;
    $p=$M,$q=$L if/O/;
    push@s,[$M,$L] if/#/;
  }
  $/=$\;
  <>
}

for(@s){
  $x=$_->[0]-$p,$y=$_->[1]-$q;  

  ($x<0&&$x<$y&&$y<=-$x ? $l :
   $x>0&&-$x<=$y&&$y<$x ? $r :
   $y<0&&$y<=$x&&$x<-$y ? $u : $d)=_;
}

$_=$l?" #":_;
$u && print;
print$l."O$r";
$d && print;